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Quickest  detection  of  the  onset  of  a  signal  is  a  common  problem  in  many  applications.  For  example, 
consider  the  detection  of  a  sonar  contact  as  it  enters  the  sonar’s  detection  range.  While  Page’s  test  is  known 
to  optimally  provide  the  lowest  average  delay  before  detection  (i.e.,  D  or  detection  latency)  for  a 
constrained  average  time  between  false  alarms,  it  does  not  necessarily  maximize  the  probability  of 
detecting  (Pd)  an  ephemeral  signal  (e.g.,  a  sonar  contact  passing  through  a  convergence  zone).  In  such 
cases  a  common  alternative  is  the  sliding  m-of-n  detector  where  a  detection  is  declared  when  m  successes 
are  observed  within  n  of  the  most  recent  trials  (e.g.,  3  detections  during  the  5  most  recent  pings). 
Techniques  for  evaluating  or  approximating  the  performance  measures  of  the  sliding  m-of-n  detector  are 
developed  and  used  to  optimally  design  the  detector.  As  expected,  Page’s  test  outperforms  the  sliding 
m-of-n  detector  with  respect  to  D  ,  except  under  certain  cases  of  signi  cant  mismatch  between  the  assumed 
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outperforms  Page’s  test  with  respect  to  Pd  or  robustness  in  false-alarm  rate  to  mismatch  in  design 
assumptions.  Unfortunately,  optimization  requires  dierent  (m;  n)  pairs  as  a  function  of  signal  length.  Thus, 
while  Page’s  test  remains  the  most  desirable  detector  to  minimize  D  or  if  the  signal  length  is  unknown,  the 
gains  in  Pd  achievable  by  a  properly  designed  sliding  m-of-n  detector  make  it  the  best  choice  for 
nite-duration  signals  of  known  length. 
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Abstract 


Quickest  detection  of  the  onset  of  a  signal  is  a  common  problem  in  many  applications.  For  example, 
consider  the  detection  of  a  sonar  contact  as  it  enters  the  sonar’s  detection  range.  While  Page’s  test 
is  known  to  optimally  provide  the  lowest  average  delay  before  detection  (i.e.,  D  or  detection  latency) 
for  a  constrained  average  time  between  false  alarms,  it  does  not  necessarily  maximize  the  probability 
of  detecting  (Pd)  an  ephemeral  signal  (e.g.,  a  sonar  contact  passing  through  a  convergence  zone). 
In  such  cases  a  common  alternative  is  the  sliding  m-of-n  detector  where  a  detection  is  declared 
when  m  successes  are  observed  within  n  of  the  most  recent  trials  (e.g.,  3  detections  during  the 
5  most  recent  pings).  Techniques  for  evaluating  or  approximating  the  performance  measures  of 
the  sliding  m-oi-n  detector  are  developed  and  used  to  optimally  design  the  detector.  As  expected, 
Page’s  test  outperforms  the  sliding  m-of-n  detector  with  respect  to  D,  except  under  certain  cases 
of  significant  mismatch  between  the  assumed  and  actual  single-trial  success  probability.  However, 
for  finite-duration  signals,  the  sliding  m-of-n  detector  outperforms  Page’s  test  with  respect  to  Pd 
or  robustness  in  false-alarm  rate  to  mismatch  in  design  assumptions.  Unfortunately,  optimization 
requires  different  (m,  n)  pairs  as  a  function  of  signal  length.  Thus,  while  Page’s  test  remains  the 
most  desirable  detector  to  minimize  D  or  if  the  signal  length  is  unknown,  the  gains  in  Pd  achievable 
by  a  properly  designed  sliding  ?n-of-n  detector  make  it  the  best  choice  for  finite-duration  signals  of 
known  length. 
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1  Introduction 


Several  detection  problems  in  sonar  signal  and  information  processing  can  be  characterized  as 
detecting  the  onset  of  a  signal  as  quickly  as  possible.  A  detection  example  is  found  in  detecting 
the  presence  of  an  active  sonar  contact  over  multiple  pings  as  it  enters  within  the  sonar’s  detection 
range.  Tracking  algorithms  provide  examples  both  in  track  initialization,  which  occurs  when  energy 
is  first  observed  consistently  along  paths  representative  of  constant-velocity  target  motion,  and  also 
when  tracks  are  ended.  Page’s  test  [1]  is  the  most  common  choice  for  the  quickest  detection  of  signal 
onset,  owing  to  its  ease  of  implementation  and  optimality  in  providing  the  minimum  average  delay 
to  detection  ( D )  for  a  constrained  average  time  between  false  alarms  (T).  A  common  alternative, 
however,  is  the  sliding  m-of-n  detector  [2,  3]  which  declares  a  detection  when  m  of  the  previous 
n  trials  are  successes  (e.g.,  a  signal  observed  on  3  of  the  5  most  recent  pings).  Noting  that  the 
optimality  of  Page’s  test  does  not  extend  to  the  probability  of  detecting  a  finite  duration  signal,  it 
is  possible  that  the  sliding  m-of-n  detector  is  a  desirable  alternative. 

The  intuitiveness  of  the  m-of-n  detection  criteria  makes  it  popular,  despite  the  limited  availability 
of  performance  assessment  techniques  when  it  is  applied  sequentially  to  detect  the  onset  of  a  signal. 
The  m-of-n  detector  is  most  often  analyzed  as  it  is  applied  in  data  fusion  [4,  Sect.  3.4];  that  is,  in  the 
context  of  m  successes  within  a  single  sample  of  n  trials.  When  sliding  the  m-of-n  detector  along 
a  data  sequence,  it  can  be  described  as  a  finite-state  Markov  process  (FSMP)  [2,3].  Williams  [3] 
developed  an  algorithm  to  compute  the  probability  of  detection  when  the  input  data  had  a  time- 
varying  single-trial  probability  of  success.  Bar-Shalom  [2,  Sect.  3.3]  presents  results  for  small  values 
of  m  and  n.  Unfortunately,  the  FSMP  characterization  encounters  computational  limitations  as  n 
increases.  The  sliding  m-of-n  detector  can  also  be  described  as  a  scan  statistic  [5].  Approximations 
for  the  average  stopping  time  of  the  test  can  be  found  in  [5,  Sect.  4.2];  however,  they  assume  a 
minimum  decision  time  of  n  samples  and  are  therefore  most  appropriate  for  use  under  false  alarm 
conditions.  The  focus  of  this  report  is  on  techniques  for  evaluating  performance  measures  for  the 
sliding  m-of-n  detector,  designing  the  detector  to  meet  a  false-alarm  rate  (FAR)  specification  and 
minimize  D  or  maximize  the  probability  of  detection  {Pd)  of  a  finite-duration  signal,  and  to  compare 
the  performance  of  the  sliding  m-of-n  detector  to  that  of  Page’s  test. 

The  m-of-n  detector  operates  on  binary  data  which  are  typically  obtained  through  an  interme¬ 
diate  thresholding  process.  In  Sect.  2,  the  quickest  detection  performance  measures  and  their 
approximations  for  Page’s  test  are  reviewed  and  applied  to  exponentially-distributed  data  and  the 
Bernoulli-distributed  binary  data  they  produce  upon  thresholding.  Analysis  of  the  quickest  de¬ 
tection  receiver  operating  characteristic  (ROC)  curves  illustrates  the  performance  lost  from  using 
thresholded  data  rather  than  pre-thresholded  data.  Approximations  are  presented  for  the  data 
threshold  maximizing  the  asymptotic  efficiency  of  Page’s  test  for  several  common  signal  and  noise 
scenarios. 
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In  Sect.  3,  the  standard  quickest  detection  performance  measures  ( D  and  T)  are  derived  for  the 
sliding  m-of-n  detector  using  Williams’  FSMP  description  and  stationary  data  statistics.  Unfor¬ 
tunately,  the  number  of  states  in  the  FSMP  characterization  of  the  sliding  m-of-n  detector  grows 
exponentially  with  n,  which  makes  analysis  difficult  for  n  >  10.  Approximations  are  therefore 
developed  that  are  applicable  to  T  and  D,  amenable  to  evaluation  for  much  larger  values  of  n  than 
the  FSMP  exact  solutions,  and  used  to  optimally  choose  the  (m,  n )  pair  meeting  a  FAR  specifica¬ 
tion  and  minimizing  D.  A  technique  is  also  presented  to  choose  the  (m,  n)  pair  meeting  the  FAR 
specification  and  approximately  maximizing  the  probability  of  detection  (P^)  of  a  finite-duration 
signal  using  an  easily  evaluated  lower  bound  on  P, 1  for  the  m-of-n  detectors. 

Finally,  the  performance  of  Page’s  test  and  the  sliding  ?n-of-n  detector  is  compared  in  Sect.  4  in 
terms  of  average  delay  to  detection  and  probability  of  detection.  Conclusions  are  presented  in  Sect. 
5  and  Matlab  code  for  some  of  the  m-of-n  subroutines  is  found  in  Appendix  B. 


2  Quickest  detection  and  Page’s  test 

Let  the  input  data  be  the  sequence  of  independent  samples 

•  •  •  -Afc— i,  Afc,  , . . .  (1) 

where  for  k  <  ko  the  data  follow  probability  density  function  (PDF)  fo{x)  and  PDF  fi(x)  for 
k  >  ko .  The  time  of  the  change  or  signal  onset,  ko ,  is  unknown.  The  quickest  detection  problem 
is  to  declare  that  the  change  has  occurred  using  the  fewest  number  of  samples  beyond  ko-  Page’s 
test  [1]  is  one  of  the  most  common  detection  algorithms  for  the  quickest  detection  of  the  onset  of 
a  signal.  Its  detection  statistic  is  obtained  through  the  simple  update 

Yk  =  max  {0,  Yk_{  +  g  (Xk)\  (2) 

with  the  initial  value  taking  on  Yq  =  0.  When  the  function  g{x )  is  the  log-likelihood  ratio  (LLR), 

g(x)  =  log  yyU  ’  (3) 

L/oOuJ 

no  other  detection  algorithm  provides  a  smaller  average  delay  before  detection  while  simultaneously 
meeting  a  constraint  on  the  average  time  between  false  alarms  [6,7]. 

For  the  purposes  of  comparing  Page’s  test  with  a  sliding  m-of-n  detector,  performance  measures 
for  quickest  detection  are  reviewed  and  evaluated  for  Page’s  test.  An  exponentially-distributed 
data  example  is  used  to  produce  the  Bernoulli-distributed  (i.e.,  binary)  data  on  which  the  m-of-n 
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detector  operates.  The  binary,  decision-level,  data  are  produced  by  the  exponentially-distributed 
data  through  an  intermediate  thresholding  process  (e.g.,  as  in  de-centralized  detection  where  only 
local  decisions  are  transmitted  to  a  fusion  agent). 


2.1  Performance  measures 


While  the  probability  of  detection  (P^)  and  probability  of  false  alarm  (P/a)  are  the  standard 
performance  measures  for  hypothesis  testing  with  fixed  sample  sizes,  the  quickest  detection  problem 
requires  different  measures.  Under  the  null  hypothesis,  the  average  time  between  false  alarms  (T) 
is  used  to  represent  how  false-alarm  prone  a  detector  is  while  the  average  delay  to  detection  (D) 
quantifies  detection  performance  under  the  alternative  hypothesis  of  signal  presence.  T  and  D  are 
sometimes  called  the  average  stopping  times  (ASTs),  average  sample  numbers  (ASNs)  or  average 
record  lengths  (ARLs).  In  terms  of  performance  measures  of  interest  to  sonar  signal  processing, 
T  is  the  inverse  of  the  false  alarm  rate  while  D  is  detection  latency.  These  measures  suffice  when 
the  duration  of  the  signal  is  infinite  once  it  starts  (e.g.,  as  in  machinery  failure)  or  at  least  much 
longer  than  D.  However,  they  can  be  supplemented  with  when  the  signal  is  ephemeral;  see  [8] 
for  several  techniques  to  evaluate  the  probability  of  detecting  finite-duration  signals  using  Page’s 
test. 


The  detector  data  mapping  function  g(x)  is  often  not  the  LLR  which  requires  some  knowledge 
about  the  signal  (e.g.,  its  strength).  When  the  LLR  is  used,  it  satisfies  the  inequalities 


E0  [g(X)}  <  0  <  Pi  [g(X)} 


(4) 


where  Po[-]  and  Pi[-]  represent  taking  the  expectation  over,  respectively,  the  null  and  alternative 
hypotheses.  Thus,  when  no  signal  is  present  Page’s  test  is  often  reset  to  zero,  but  then  tends  to  rise 
when  signal  is  present.  When  g(x )  is  not  the  LLR,  good  performance,  as  defined  by  an  exponential 
relationship  between  the  threshold  h  (to  which  Yj.  is  compared  in  declaring  a  detection)  and  T  and 
a  linear  relationship  between  h  and  D,  requires  the  inequalities  in  (4)  be  satisfied.  When  they  are, 
approximations  to  the  ASNs  are  [9,  Ch.  5] 


hto 


and 


D 


1  +  hto  ~  e 
t0E0  [g{X)\ 

1  +  ht\  -  ehtl 


(5) 


(6) 


tiEi  [g(X)\ 

where  to  and  t\  are  the  non-zero  unity  roots  of  the  moment  generating  function  (MGF)  of  g(X) 
That  is,  U  /  0  such  that 


Ei 


Ji9(X) 


=  1 


(7) 
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for  i  =  0  and  1. 

In  most  practical  situations,  T  will  be  large,  which  will  require  a  large  threshold  h.  In  this  regime, 
and  under  the  conditions  described  in  (4),  the  asymptotic  efficiency  {if)  of  Page’s  test  [10]  describes 
the  log-linear  relationship  between  T  and  D , 


E-s 

■e 

II 

c- 

(8) 

=  toEMX)]. 

(9) 

2.2  Exponentially-  and  Bernoulli-distributed  data  examples 

The  primary  focus  of  this  paper  is  on  quickest  detection  of  Bernoulli-distributed  data.  However, 
it  is  important  to  recall  the  inherent  loss  in  performance  incurred  by  working  with  such  a  harsh 
quantization.  To  illustrate  the  loss  in  performance,  an  exponentially-distributed-data  example  is 
used  to  derive  the  Bernoulli-distributed  data  through  the  standard  thresholding  process.  In  order 
to  evaluate  performance,  the  data  mapping  function  g{x )  must  be  derived  along  with  the  MGF 
unity  roots  and  means  under  the  null  and  alternative  hypotheses.  While  g{x )  and  its  mean  values 
are  evaluated  here,  the  procedure  for  obtaining  the  MGF  unity  roots  is  described  in  Appendix  A. 


2.2.1  Exponentially-distributed  data 

Suppose  the  data  follow  an  exponential  distribution  with  mean  A  under  the  null  hypothesis  and 
mean  A(1  +  s)  under  the  alternative.  Thus,  A  is  the  average  noise  power  and  s  is  the  signal-to-noise 
power  ratio  (SNR).  The  LLR  is 

Kx)  =  vyy-  [t  -  (l+  log (1  +  s)l  (10) 

1  +  s  [A  V  s J 

and  clearly  depends  on  both  s  and  A.  In  this  paper  A  is  assumed  to  be  known,  although  it  is 
typically  estimated  in  a  normalization  procedure.  The  detector  is  implemented  with  a  design  SNR, 
s,  which  will  most  likely  differ  from  the  observed  SNR  s.  The  optimality  of  the  LLR  in  Page’s  test 
is  not  altered  by  a  change  in  scale,  so  the  data  mapping  function 

g(X)  b(s)  (li) 

is  used  where  the  bias 

b{$)  =  ^1  +  l°g(l  +  s)  (12) 
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only  depends  on  the  design  SNR. 


The  mean  of  the  data  mapping  function  is 


under  the  null  hypothesis  and 


under  the  alternative  hypothesis. 


E0[g(X)]  =  l-b(s) 


Ei[g(X)\  =  1  +  s-b(s) 


(13) 

(14) 


2.2.2  Bernoulli-distributed  data 


If  the  exponentially-distributed  data,  normalized  by  the  noise  power,  are  compared  with  a  threshold 
( hx )  prior  to  application  to  Page’s  test,  the  binary  data  sequence  can  be  described  as 

=  (15) 

where  /(•)  is  an  indicator  function  returning  one  when  the  argument  is  true  and  zero  otherwise. 
The  datum  is  a  Bernoulli  random  variable  with  probability 

Po  = 


under  the  null  hypothesis  (Hq)  and 

Pi  = 


under  the  alternative  hypothesis  {Hi).  Note  that  these  are  the  well  known  false-alarm  and  detection 
probabilities  for  envelope  detection  of  a  Gaussian  target  in  Gaussian  noise. 

Scaling  the  LLR  of  the  Bernoulli  data  to  isolate  the  dependence  on  po  and  p\  into  the  bias  results 
in  the  data  mapping  function 

g{u)  =  u  -  b{p0,pi)  (18) 

where  the  bias  is 

(19) 
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b(po,Pi )  = 


1  + 


log(f>i/po) 


log[(l  Po )  /  ( 1  -  .Pi)] 


Pr  >  hx  \H0 

^  hx 


Pr  <!  >  hx  |  Hi 


e  i+s  =p0 


1 

1+s 


(16) 


(17) 
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If  the  detection  process  producing  the  Bernoulli  data  is  unbiased,  then  p\  >  po  and  b  ( po,Pi )  E  (0, 1). 
While  it  is  reasonable  to  expect  po  to  be  known,  it  is  unlikely  that  p\  will  be  known  with  any 
accuracy.  As  was  done  for  the  exponentially-distributed  data,  design  probabilities  po  and  pi  can 
be  used  to  form  the  bias 


b  =  b(po,pi) 


(20) 


or  b  may  be  chosen  in  (0, 1). 


The  mean  value  of  the  data  mapping  function  is 

Eo  [ g(U )]  =  P0  -  b 

(21) 

under  the  null  hypothesis  and 

Ei  [g(U)\  =pi-b 

(22) 

under  the  alternative  hypothesis. 

2.2.3  Performance  comparison 


The  receiver  operating  characteristic  (ROC)  curves  for  Page’s  test  operating  on  the  exponentially- 
and  Bernoulli-distributed  data  are  compared  in  Fig.  1  for  the  case  of  po  =  10-4  and  p\  =  0.5, 
which  equates  to  s  ~  10.9  dB  for  the  exponentially-distributed  data  model.  The  solid  lines  in  the 
figure  are  D  as  a  function  of  log10  T  when  the  bias  is  chosen  correctly  according  to  s  or  p\.  The 
asymptotic  efficiencies  (as  computed  via  (9))  for  these  cases  are 

tjb  =  3.9  and  tje  =  9.7  (23) 

for  the  Bernoulli-  and  exponentially-distributed  data,  respectively.  Their  ratio  tje/ob  =  2.5  is 
representative  of  the  factor  by  which  D  for  the  exponentially-distributed  data  is  below  that  for  the 
Bernoulli-distributed  data;  that  is,  the  signal  is  detected  approximately  2.5  times  faster  using  the 
original  data  compared  to  using  the  quantized,  binary  data. 

The  dashed  lines  in  Fig.  1  represent  the  performance  when  the  bias  in  the  data  mapping  function 
is  chosen  incorrectly.  For  the  exponential  data,  the  bias  was  chosen  assuming  an  SNR  of  15,  5,  2, 
and  0  dB,  with  performance  increasingly  worse  in  this  order.  Although  not  always  the  case,  this 
illustrates  that  overestimating  the  SNR  is  often  not  as  bad  as  overestimating  it.  The  dashed  curves 
for  the  Bernoulli  data  are  formed  with  a  design  p±  formed  from  the  design  SNR  according  to 

Pi  =  P(P  (24) 

and  illustrate  similar  performance  degradation  relative  to  the  correctly  chosen  bias,  which  optimizes 
performance. 
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Figure  1:  ROC  curves  for  Page’s  test  operating  on  exponential  and  Bernoulli  data  with  correct  bias 
(solid  lines)  and  a  mismatched  bias  (dashed  lines).  The  advantage  of  working  with  data  prior  to 
quantization  is  evident  in  the  lower  values  of  D  for  the  exponential  data. 
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2.3  Optimal  quantization 


When  the  pre-thresholded  data  are  accessible  or  the  quantization  process  is  included  in  the  detector 
design,  the  quantization  threshold  hx  can  be  chosen  to  maximize  the  asymptotic  efficiency  of  Page’s 
test  formed  from  the  binary  data  of  (15).  When  the  FAR  specification  results  in  a  high  threshold 
in  Page’s  test,  then  optimizing  i]  will  also  minimize  D  for  a  given  T. 

Using  the  definition  of  p  in  (9),  the  optimal  bias  from  (19),  and  the  approximation  for  to  f°r 
Bernoulli-distributed  data  found  in  Appendix  A, 


yields  the  approximation 


to 


-  log  (go) 

1-6 


n 


-Pi  log  (po)  + 


(1  -pi)log(j=2a) 

log(Pl)  _  ^ 
log(po) 


(25) 


(26) 


For  the  fluctuating  signal  described  in  Sect.  2.2.1,  the  first  term  in  (26)  dominates  the  second, 
especially  as  SNR  increases  (where  p\  — >  1).  Using  (16)  and  (17)  for  po  and  pi,  the  first  term  in 
(26)  is  seen  to  be  maximized  when 

hx  =  s  -\-  1.  (27) 

Starting  with  this  threshold  and  implementing  one  step  of  a  Newton-Raphson  update  to  find  the 
maximum  of  (26)  results  in  the  refined  approximation  to  the  optimal  threshold  (in  dB), 


HI 


10  log10  hi 


10  log10(s  +  1)  +  10  log10 


1  + 


s-(e-  l)-1 
1  +  log  (1  -  1/e) 


-i' 


(28) 


where  e  =  e1  «  2.7183. 


The  approximation  of  (28)  is  compared  in  Fig.  2  with  the  optimum  threshold  for  the  fluctuating 
signal  as  obtained  through  a  numerical  search  using  (9).  As  seen  in  the  figure,  the  approximation 
is  quite  good,  especially  for  larger  values  of  SNR.  The  average  absolute  error  is  only  0.026  dB  for 
SNR  ranging  from  2  to  16  dB.  Although  the  maximum  error  is  more  than  one-tenth  of  a  dB,  it 
occurs  for  the  smallest  SNR.  Three  other  common  signal  and  noise  scenarios  are  also  considered:  a 
non-fluctuating  signal  in  complex  Gaussian  noise,  a  Gaussian  shift-in-mean  signal,  and  a  Gaussian 
change-in- variance  signal.  The  statistical  characterizations  of  the  detection  statistic  under  the  null 
and  alternative  hypotheses  along  with  equations  for  po  and  p\  for  the  various  scenarios  considered 
are  found  in  Table  1. 
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Table  1:  Data  statistical  characterizations  under  the  null  and  alternative  hypotheses  and  the  asso¬ 
ciated  success  probabilities  for  various  signal  and  noise  scenarios  as  a  function  of  the  data  threshold 
hx  and  SNR.  An  exponential  distribution  with  mean  A  is  denoted  by  £(A),  a  non-central  chi-squared 
distribution  with  two  degrees  of  freedom  and  non-centrality  parameter  5  by  xi(^)>  an(A  a  Gaussian 
distribution  with  mean  p  and  variance  a2  by  Af(p,  a2).  The  cumulative  distribution  functions 
(CDFs)  of  the  non-central  chi-squared  and  Gaussian  distributions  are,  respectively,  Fx 2  (x',5)  and 
3 >{(x-  n)/cr). 


Signal  type 

H0 

Po 

Pi 

Fluctuating 

X~£(l) 

e~hx 

X  ~  £(1  +  s ) 

e~hx/(l+s) 

Non-fluctuating 

X~£(l) 

e~hx 

2 A  ~  x!(2 a) 

1  -  FX2  (2 hx;  2s) 

Gaussian  shift-in-mean 

x  ~  jV(o,  1) 

1  -  (hx) 

x  ~  N(i/s,  1) 

1  -$(hx-  v/s) 

Gaussian  change-in-variance 

X~Af(0,l) 

1  -  <I>  (. hx ) 

X  ~  AA(0, 1  +  s) 

1  —  <h  (hx/yj  1  +  s) 

The  complicated  forms  of  the  threshold/exceedance-probability  relationships  for  these  latter  signal 
types  hinders  derivation  of  an  approximation  similar  to  (28).  However,  noting  the  smooth  functional 
forms  of  the  optimal  threshold  seen  in  Fig.  2  on  logarithmic  scales,  the  functions  were  fit  via  least- 
squared-error  optimization  using  half-powers  of  SdB  =  10 logit)5  (i-e->  V'S'dB,  •SdB,  ancl  'S’dEs)-  The 
functional  forms  of  the  approximations  to  the  optimal  threshold  are  shown  in  Table  2  for  each 
signal  and  noise  scenario  along  with  the  average,  median  and  maximum  errors  over  the  range  of 
SNR  evaluated.  A  least-squared-error  fit  to  the  optimal  threshold  for  the  fluctuating  signal  (not 
shown  in  Fig.  2)  provided  a  better  overall  fit  than  (28)  at  the  expense  of  an  increase  in  error  for 
large  SNR.  While  the  Gaussian  signals  only  required  two  terms  to  keep  the  maximum  absolute 
error  less  than  one-twentieth  of  a  decibel  over  the  range  of  SNR  evaluated,  the  fluctuating  and 
non-fluctuating  signals  required  four  terms.  As  seen  in  Fig.  2,  the  fit  is  quite  good  within  the  range 
of  evaluation.  Other  than  (28),  these  functions  are  not  expected  to  work  beyond  S^b  £  [2, 16]. 

The  general  trend  for  the  optimal  threshold  is  one  that  increases  with  SNR,  indicating  that  7/  is 
maximized  primarily  by  reducing  po,  although  p\  was  also  slightly  increasing  with  SNR.  The  errors 
in  the  approximations  found  in  Table  2  are  small  enough  that  they  can  be  used  to  design  the 
intermediate  thresholding  when  it  is  desired  to  apply  Page’s  test  to  the  binary,  decision-level  data. 
As  will  be  seen  in  the  following  section,  this  is  likely  a  good  choice  when  an  m-oi-n  detector  is  used 
instead  of  Page’s  test. 
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SNR  (dB) 


Figure  2:  Data  threshold  (hx)  optimizing  the  asymptotic  efficiency  of  Page’s  test  for  various  signal 
and  noise  scenarios  (solid  fines).  The  general  trend  observed  is  for  an  increasing  threshold  with 
SNR.  Approximations  to  the  optimal  threshold  are  shown  as  dashed  lines  or  ‘+’  marks. 
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3  Sliding  m-of-n  detectors 


A  common  technique  for  fusing  information  is  to  apply  an  m-of-n  rule  where  a  detection  declaration 
is  made  whenever  m  successes  occur  within  n  associated  trials.  The  situation  is  common  in  dis¬ 
tributed  detection  [4,  Sect.  3.4]  and  known  as  binary  integration  in  the  radar  community  [11,  Sect. 
6.4]  owing  to  its  equivalence  to  thresholding  a  rolling  sum  of  binary  data  representing  the  individual 
trials.  Analysis  of  the  m-of-n  detector  is  usually  limited  to  the  static  case  where  only  one  set  of 
n  inputs  is  evaluated.  Sliding  an  m-of-n  detector  along  a  data  sequence,  which  is  often  how  it  is 
used  to  reduce  false  alarms,  results  in  what  is  known  in  the  statistics  literature  as  a  scan  statis¬ 
tic  [5].  It  can  also  be  described  as  a  discrete-time,  finite-state  Markov  process  (FSMP)  [2,3].  In  [3] 
a  technique  was  presented  for  numerically  evaluating  the  probability  of  detection  for  an  input  se¬ 
quence  of  Bernoulli  random  variables  with  time- varying  probabilities — the  example  considered  in  [3] 
was  an  active  sonar  contact  transiting  through  an  area  presenting  a  limited  detection  opportunity. 
BarShalom  [2,  Sect.  3.3]  presents  results  for  small  values  of  m  and  n. 

In  this  section,  the  quickest-detection  performance  measures  ( D  and  T )  are  developed  along  with 
a  more  simply  evaluated  alternative  to  Williams’  algorithm  for  Pd  and  a  lower  bound  on  Pcj  for  the 
stationary-data  case  in  preparation  for  comparing  the  performance  of  the  sliding  m-of-n  detector 
with  Page’s  test.  While  the  characterization  of  the  sliding  m-of-n  detector  as  an  FSMP  allows 
for  exact  evaluation  of  the  average  stopping  times  and  Pd,  it  is  limited  by  exponential  growth  in 
the  number  of  states  as  n  increases.  Approximations  to  the  average  stopping  time  of  the  scan 
statistic  [5,  Sect.  4.2]  surmount  this  limitation,  but  assume  a  minimum  decision  time  of  n  rather 
than  m  and  have  numerical  evaluation  issues  for  small  p.  To  surmount  this  limitation,  a  small-p 
approximation  to  T  is  presented  along  with  an  alternative  computational  method  approximating 
the  average  stopping  times  that  is  useful  for  a  much  wider  range  of  n.  Finally,  techniques  are 
presented  for  choosing  m  and  n  to  meet  a  false  alarm  specification  and  to  approximately  maximize 
the  probability  of  detecting  a  finite  duration  signal  and  minimizing  the  average  delay  to  detection. 


3.1  m-of-n  detection  as  a  finite-state  Markov  process 

The  key  to  Williams’  analysis  in  [3]  was  to  describe  the  m-of-n  detector  as  an  FSMP  with  the 
states  defined  by  the  previous  n  —  1  trials,  what  he  termed  the  detector’s  history.  The  detection 
event,  observing  m  successes  in  n  consecutive  trials,  was  characterized  by  an  ‘accepting’  state  while 
all  other  states  were  defined  as  non-accepting.  The  complete  set  of  states  was  enumerated  by  the 
integer  value  of  the  binary  number  defined  by  the  history.  Consider,  for  example,  m  =  3  and 
n  =  5.  If  the  previous  4  trials  were  the  sequence  failure,  success,  success,  failure  then  the  history 
is  the  binary  number  0110  with  a  state  enumeration  of  the  integer  6.  If  the  next  trial  is  a  success, 


UNCLASSIFIED 


approved  for  public  release 


12 


CausaSci  LLC 


UNCLASSIFIED 


2011-01 


then  this  state  transitions  to  the  accepting  state  having  observed  3  success  in  5  consecutive  trials. 
However,  a  failure  in  the  next  trial  would  transition  to  state  1100  (integer  value  12).  A  technique 
for  obtaining  an  equivalent  FSMP  with  a  reduced  set  of  states  was  also  presented  in  [3],  enabling 
a  more  efficient  numerical  analysis  for  larger  values  of  n.  Matlab  code  for  creating  the  structures 
defining  the  Markov  process  is  presented  in  Appendix  B. 

The  Markovian  nature  of  the  process  implies  that  the  vector  p(k)  of  probabilities  of  being  in  each 
of  the  states  at  time  k  can  be  determined  from  the  probabilities  at  time  k  —  1  according  to  [12,  Sect. 
4.3] 

pT(k )  =  pT(k  —  l)P(k)  (29) 

where  P(k )  is  the  (potentially)  time- varying  state  transition  matrix  and  the  superscript  T  represents 
the  transpose  operation.  Williams  presented  in  [3]  an  algorithmic  update  accomplishing  (29); 
however,  the  matrix  representation  is  maintained  here  owing  to  its  simpler  implementation  and 
the  corresponding  ability  to  evaluate  D  and  T  with  a  simple  summation  over  the  FSMP  states  for 
stationary  data. 

Suppose  there  are  ns  states,  including  the  accepting  or  absorbing  state.  In  ordering  the  states  within 
p(k),  let  the  first  state  be  a  history  comprising  all  zeros  (i.e.,  none  of  the  past  n  —  1  trials  were 
successful)  and  let  the  last  state  be  the  accepting  state  (i.e.,  the  m-of-n  success  state).  Assuming 
that  once  the  accepting  state  is  entered,  the  process  remains  there,  the  state  transition  matrix  can 
be  described  by  the  decomposition 

p/m  _  [  pcc{k)  pcs{k ) 

p{k)  -  y  qT  1 

where  Pcc{k)  is  the  transition  matrix  for  the  ns  —  1  continuing  states,  pcs{k)  is  a  vector  of  the 
transition  probabilities  from  the  continuing  states  to  the  accepting  state,  and  0  is  a  vector  of  zeros. 

The  probability  of  being  in  one  of  the  continuing  states  in  the  kth  trial  can  then  be  described  by 
the  recursion 

Pc  ( k )  =  Pi ( k  -  l)Pcc(k)  (31) 

where  pc( 0)  is  a  vector  of  the  initial  probabilities  of  being  in  a  continuing  state. 

The  performance  of  the  sliding  m-of-n  detector  is  described  by  defining  the  probability  mass  func¬ 
tion  (PMF)  of  the  stopping  time  K ,  which  is  defined  as  the  first  time  the  process  enters  the 
accepting  state.  The  PMF  for  K  is  the  probability  of  transitioning  to  the  accepting  state  at  time 
k  from  one  of  the  continuing  states  at  time  k  —  1, 

fK{k)  =p^{k-  l)pcs(k)  (32) 
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for  k  =  1,  2,  . . . . 

The  probability  of  detecting  a  signal  within  k  trials  using  a  sliding  m-of-n  processor  is  the  CDF  of 

K, 

k 

Pd{k)  =  FK{k)  =  YJfK(l).  (33) 

i=i 

This  also  approximates  the  probability  of  detecting  a  finite-duration  signal  of  length  k,  only  ignoring 
any  detections  occurring  after  time  k  but  before  time  k  +  n  (i.e.,  while  the  sliding  window  is  still 
affected  by  the  signal).  For  Page’s  test,  these  were  termed  latent  detections  [8]. 

The  average  sample  numbers  (T  and  D)  are  the  mean  of  K  under  the  null  and  alternative  hypothe¬ 
ses,  respectively.  These  are  typically  evaluated  when  the  process  is  stationary  (i.e.,  pk  is  constant) 
and  therefore  discussed  after  the  simplifications  of  the  following  section. 

The  time-varying  part  of  the  transition  matrix  and  vector  is  limited  to  the  success  probability  of 
the  kth  trial.  Defining  this  as  pk,  the  continuing  state  transition  matrix  can  be  described  as 

Pcc(k)  =PkJ+  +  (1  ~Pk)J-  (34) 

where  J  +  and  J_  are  indicator  matrices  (the  elements  containing  either  a  zero  or  a  one)  representing 
which  state  is  achieved  at  the  next  time  step  given  a  success  («7+)  or  a  failure  («/_).  Similarly,  the 
probability  of  transitioning  from  the  continuing  to  accepting  states  can  be  described  by 

Pcs(k)=Pkj+  (35) 

where  j+  is  an  indicator  vector.  Using  (34)  and  (35)  in  (31)  and  (32)  is  more  efficient  in  Matlab  than 
the  algorithm  developed  in  [3].  Matlab  code  for  producing  the  indicator  matrices  from  Williams’ 
Markov-process  structure  definition  is  found  in  Appendix  B. 


3.2  State  and  stopping  time  probabilities  for  stationary  data 

When  the  state  transition  probabilities  do  not  change  with  time,  simplifications  allow  further 
analysis  of  the  m-of-n  processor.  Dropping  the  dependence  on  k  of  Pcc,  pcs,  and  p,  the  update 
equation  (31)  for  the  continuing  state  probabilities  becomes 

PTA.k)  =  pJ(0)Pkc  (36) 

=  pJ{0)VAkV~1 

where  V  and  A  are  the  eigenvector  and  eigenvalue  matrices  of  Pcc.  Noting  that  A  is  diagonal 
illustrates  the  ease  with  which  the  state  probabilities  may  be  calculated. 
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Similarly,  the  PMF  for  K  is  simplified  to 

fK(k)  =  pno)-PcC_ipcs 

=  pT(0)V  A^V-'p^ 

ns  — 1 

=  £  oAXt1  (37) 

i= 1 

where  A*  are  the  eigenvalues  of  Pcc  (i.e.,  the  diagonal  elements  of  A)  and  a*  and  bi  are  elements  of 
the  vectors 

a  =  Urpc(0)  (38) 

and 

b  =  V~1Pcs.  (39) 

The  CDF  of  K  (yielding  Pd(k))  can  also  be  evaluated  without  resorting  to  the  cumulative  summa¬ 
tion  of  (33), 


FK(k) 


k 


y  /k(i) 


na— 1  k 


i=  1  1=1 


na— 1 

y  ciA 

i=  1 


(!  -  A*) 
(1  ~  A*)  ' 


(40) 


Note  that  the  PMF  of  the  stopping  time  for  times  k  =  m,. . . ,  n  is  simply  the  probability  of 
observing  m  —  1  successes  in  k  —  1  trials  followed  by  a  success, 


/a'(^) 


k-  1 
m  —  1 


Pm(l-P) 


k—m 


(41) 


where  ‘i  choose  j\ 

(i\_  *!  _  r(*  +  l) 

W  “  j!(7  -  j)!  "  P(j  +  l)r(i  -  i  +  1)  ’ 

is  the  number  of  different  combinations  of  size  j  from  a  population  of  size  i. 
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The  CDF  of  the  stopping  time  is  simply  the  summation 

=  £  („V-1i)p’”(1  -  p)>~m  (43) 

j =m  '  / 

and  provides  an  easily  evaluated  alternative  to  the  FSMP  result  of  (40)  that  is  extremely  useful 
when  n  is  large. 


The  initial  state  probability  vector  for  the  continuing  states  (pc(0))  should  be  chosen  assuming 
the  sliding  m-of-n  detector  has  been  operating  for  some  time  under  the  null  hypothesis,  but  not 
resulted  in  a  false  alarm.  Following  [8],  the  update  equation  for  the  continuing  state  probabilities 
given  no  false  alarms  is 


Pc\dk  ~  l)pcc 
PTc\c{k  -  1)-PCC1 


(44) 


where  1  is  a  vector  of  ones.  Using  the  eigen-decomposition  of  (36),  it  is  straightforward  to  show 
that  the  limiting  continuing  state  probability  vector  is  the  left  eigenvector  of  Pcc  associated  with 
its  largest  eigenvalue 


eiV 


-i 


zTV~ll 


(45) 


where  e\  is  a  unit  vector  pointing  in  the  direction  of  the  first  element  and  it  is  assumed  that  the 
largest  eigenvalue  is  in  the  first  position. 


For  small  values  of  p,  an  accurate  approximation  can  be  obtained  from  the  probability  of  observing 
the  number  of  successes  in  the  history  of  each  continuing  state.  Defining  (3(i )  as  the  number  of 
successes  in  the  length  n  —  1  history  of  the  ith  state,  the  asymptotic  probability  of  observing  the 
ith  continuing  state  is  approximately 


Pc|cM 


ns~l  -nPti) 


Ens- 
7=1 


ph 


[  pPW  pP(2)  ...  pO(ns-l)  jT 


(46) 


For  values  of  n  <  10,  the  approximation  of  (46)  is  quite  good  when  p  <  0.1.  Define  the  error  as 
the  maximum  over  all  values  of  m  for  a  given  n  of  the  sum  of  the  total  absolute  error  between  the 
FSMP  exact  solution  of  (45)  and  the  approximation  of  (46).  The  error  was  worst  for  small  values 
of  m  >  1  and  larger  values  of  n,  and  ranged  from  p 2  for  the  smaller  values  of  n  to  p  for  the  larger 
values.  Noting  that  the  first  state  ( i  =  1)  is  the  only  state  with  no  successes  in  the  history  (i.e., 
(3{  1)  =  0),  most  of  the  weight  is  placed  on  this  state  when  p  is  small.  Therefore,  using 


Pc\c(oo)  ~  ei  (47) 

is  a  reasonable  approximation  when  p  is  sufficiently  small. 
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3.3  Average  stopping  time 


Under  the  stationary  conditions  described  in  the  previous  section,  the  average  stopping  time  of  the 
sliding  m-of-n  detector  is 


E[K\ 


^2kfK(k) 

k=  1 


ns  —  1  oo 

E  EfcA"_1 


i= 1  fc=  1 


i= 1 


CLibi 

(1  -  A,:)2 


Pc(°)  ( I  -Pec)-2  Pcs 


(48) 


(49) 

(50) 


The  average  time  between  false  alarms  (T)  is  obtained  by  using  po  in  forming  the  state  transition 
probabilities,  while  D  requires  use  of  p\ . 


When  m  =  1,  the  stopping  time  of  the  sliding  m-of-n  detector  is  simply  the  time  of  the  first  success 
on  an  individual  trial,  which  is  a  geometric  random  variable.  The  PMF  of  the  stopping  time, 
regardless  of  the  value  of  n,  is  then 

fK(k)=p(l-p)k-1  (51) 

with  the  associated  CDF 

FK{k)  =  l-(l-p)k  (52) 

for  k  =  1,  . . . .  The  average  stopping  time  is 

E[K\  =  -  (53) 

P 

and  standard  deviation 

Std[A]  =  n/1~P.  (54) 

p 


3.3.1  Small-p  approximation 


When  the  probability  of  success  on  an  individual  trial  is  sufficiently  small,  the  sliding  m-of-n 
detector  can  be  approximated  as  a  sequence  of  independent  Bernoulli  trials  where  the  success 
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probability  is  the  probability  of  observing  a  failure,  followed  by  to  —  1  successes  in  the  next  n  —  1 
trials,  followed  by  one  success, 


V 


(1  ~P)  • 


pm-\l-p)n-m 


(n 


-  l)lpm{l-p)n-m+1 
( m  —  l)!(n  —  m)\ 


■P 


(55) 


Note  that  (55)  is  slightly  more  accurate  than  using  the  probability  of  observing  m  successes  in 
n  trials,  which  does  not  require  a  success  on  the  final  trial  or  a  failure  in  the  trial  immediately 
preceding  the  set  of  n  trials  leading  to  a  detection.  Similar  to  the  m  =  1  case,  the  stopping  time  is 
then  approximately  a  geometric  random  variable  with  statistics  characterized  by  (51)— (54)  using  p 
from  (55)  in  place  of  p.  The  average  stopping  time  is  approximated  as 


E[K\ 


(m—  1  )!(n  —  m)\ 

(n-  l)!pm(l  -p)n~m+ 1’ 


(56) 


which  is  significantly  easier  to  evaluate  and  less  prone  to  computational  errors  than  the  exact  FSMP 
solution  of  (50).  Support  for  the  geometric-distribution  assumption  comes  from  the  form  of  the 
stopping  time  PMF  seen  in  (37),  which  is  the  weighted  sum  of  the  eigenvalues  of  Pcc  raised  to  the 
power  k  —  1.  As  k  increases,  the  PMF  increasingly  becomes  dominated  by  the  contributions  of  the 
maximum  eigenvalue, 

fK(k)  «  aiM^1  (57) 

where  it  is  assumed  that  Ai  is  the  largest  eigenvalue.  When  p  is  small,  the  largest  eigenvalue  of 
Pcc  is  significantly  larger  than  the  next  largest  one,  additionally  contributing  to  the  accuracy  of 
this  approximation  as  seen  in  the  following  analysis. 


The  small-p  approximation  of  (56)  is  compared  with  the  FSMP  exact  solution  of  (50)  in  Figs.  3 
and  4  for  n  =  7  and  all  possible  values  of  m.  Note  that  the  exact  result  of  (53)  as  a  function  of  p 
is  shown  for  m  =  1.  The  FSMP  exact  solution  is  shown  as  asterisks  and  encounters  computational 
errors  when  p  is  too  small.  The  approximations  from  [5]  had  numerical  evaluation  issues  in  the 
same  places  as  the  FSMP  exact  solution  and  converge  to  n  as  o  — >  1  rather  than  m.  Figure  4 
expands  the  region  p  £  [0.01, 1]  where  the  poor  fit  of  the  small-p  approximation  as  p  increases  is 
revealed.  Defining  the  relative  error  as 

I  K-Ka 


where  K  is  the  average  stopping  time  from  (50)  and  Ka  is  the  small-p  approximation  of  (56),  the 
latter  (dashed  lines  in  Figs.  3  and  4)  is  quite  accurate  for  p  less  than  about  0.1  for  n  =  7  with  e# 
below  0.1  (10%  error)  and  decreasing  proportionately  with  p  (i.e.,  1%  error  for  p  <  10-2).  Note 
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that  the  average  stopping  time  is  shown  on  a  logarithmic  scale  in  the  figures.  The  relative  error 
observed  or  calculated  on  a  logarithmic  scale  is  less  than  that  computed  via  (58)  and  reported 
herein.  Although  not  shown,  the  regions  for  which  the  relative  error  is  below  10%  for  n  =  3  are 
p  <  0.25,  p  <  0.18  for  n  =  5,  and  p  <  0.06  for  n  =  10.  This  is  representative  of  a  constant  np 
effect  where  for  larger  n  the  approximation  requires  p  to  be  commensurately  smaller  for  accuracy. 
The  limited  results  presented  here  indicate  the  small-p  approximation  has  maximum  errors  around 
10%  when  np  <  0.7  or  1%  errors  when  np  <  0.07.  As  can  be  seen  in  Fig.  4,  the  errors  are  worst 
for  m  =  2  and  decrease  in  severity  as  m  increases  (e.g.,  the  error  is  less  than  10%  for  p  <  0.7  when 
m  =  n  =  7). 


3.3.2  Lower  bound 


Recall  the  stopping  time  PMF  from  (41)  for  times  k  <  n.  If  the  average  stopping  time  is  significantly 
less  than  n,  then  (41)  might  yield  an  accurate  approximation.  More  specifically,  by  letting  n  — >  oo, 
the  stopping  time  PMF  of  (41)  represents  a  lower-bound  on  the  stopping  time.  That  is,  the  test 
stopping  after  the  first  m  successes  will  stop  as  soon  or  sooner  than  the  test  stopping  the  first  time 
m  successes  are  observed  in  n  consecutive  trials.  Using  (41)  in  (48),  the  lower  bound  is  seen  to  be 

E[K]  > 


where  the  infinite  summation  in  (59) 

[13>  TaWe 

The  bound  is  compared  with  the  exact  FSMP  solution  in  Fig.  5.  As  might  be  expected,  the  bound 
is  tightest  for  m  <C  n  and  as  p  -*  1,  but  otherwise  a  very  loose  bound.  While  it  is  better  than  the 
small-p  approximation  near  p  =  1,  its  primary  utility  is  as  a  lower  bound  rather  than  an  accurate 
approximation  to  the  average  stopping  time. 


k=m 


k  —  1 
m  —  1 


Pm(l-P ) 


k—m 


oo 


mpm  r(j'  +  m  +  !)  (!  -  P)J 


pj  r(m  +  l)  j\ 


m 

P 


(59) 

(60) 


is  recognized  as  a  hypergeometric  function  simplifying  to 
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Numerical  problems  in 
FSMP  exact  solution 
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p  =  single-trial  probability  of  success 


Figure  3:  Comparison  of  the  FSMP  exact  solution  to  the  average  stopping  time  to  small-p  and 
alternative  computational  approximations  for  n  =  7  and  all  possible  values  of  m.  The  FSMP 
exact  solution  encounters  numerical  issues  when  p  is  small  enough  while  the  small-p  approximation 
degrades  as  p  increases. 
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p  =  single-trial  probability  of  success 


Figure  4:  Expanded  view  of  Fig.  3  in  the  region  p  E  [0.01, 1]  illustrating  the  failure  of  the  small-p 
approximation  and  the  accuracy  of  the  alternative  computational  approximation  as  p  increases. 
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Figure  5:  Average  stopping  time  as  a  function  of  p  for  n  =  10  comparing  the  exact  solution,  the 
alternative  computational  approximation,  and  the  lower  bound  for  various  values  of  m. 
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3.3.3  Alternative  computational  approximation 


As  seen  in  Fig.  3,  the  FSMP  exact  solution  of  (50)  encounters  numerical  issues  when  p  is  small.  It 
can  also  be  computationally  intensive  when  n  is  large  owing  to  the  substantial  number  of  states, 
which  for  m  near  n/2  grows  exponentially  with  n.  For  example,  there  are  up  to  253  states  when 
n  =  10,  up  to  6,436  when  n  =  15  and  184,757  when  n  =  20.  The  approximations  from  the  scan- 
statistics  literature  [5]  similarly  have  problems  for  small  values  of  p  and  then  are  derived  to  converge 
on  n  rather  than  m  as  p  — >  1,  so  are  primarily  useful  for  evaluating  T  when  p  is  not  too  small. 
An  alternative  computational  method  approximating  the  average  stopping  time  is  presented  in  this 
section  avoiding  the  numerical  issues  of  the  scan-statistic  approximations  and  the  FSMP  exact 
solution  both  in  defining  the  states  and  evaluating  E[K\.  Extensive  evaluation  of  (50)  revealed 
that  the  average  stopping  time  could  be  described  approximately  as 

Ep[K\  »  ^e(n»-l)G(Plc(m,n))  (61) 

where  the  function 

/•-log  p 

G(p,c)  =  /  (1  —  e~u)c  du  (62) 

Jo 

and  c(m,n )  is  a  term  that  requires  numerical  evaluation.  Noting  that  G(l,c)  =  0,  it  can  be  seen 
that  the  average  stopping  time  in  (61)  tends  to  m  as  p  — >•  1.  It  was  also  observed  that  c(m,n) 
increases  with  n,  but  decreases  with  m,  and  that  c(2,  n)  ~  n  —  1  for  large  n.  While  there  does  not 
appear  to  be  any  simple  closed-form  solutions  to  the  integral  in  (62),  the  argument 

9c(u)  =  (1  -  e~u)c  (63) 

is  suitable  for  numerical  integration.  The  procedure  implemented  in  the  Matlab  code  found  in 
Appendix  B  exploits  a  second-order  Taylor-series  expansion  of  gc(u)  about  uq  e  (0,  —  logp) 

9c(u)  ~ 


for  a  set  of  points  equally  spaced  within  the  integration  limits.  The  integral  in  (62)  is  then  approx¬ 
imated  as  the  sum  of  the  integrals  over  the  piece-wise  quadratic  functions. 

The  value  of  c(m,n)  is  chosen  to  minimize  the  error  between  the  numerical  evaluation  of  (61)  and 
the  small-p  approximation  of  (56)  for  some  sufficiently  small  value  of  p.  Although  the  minimization 
is  not  a  computationally  burdensome  evaluation,  c(m,  n)  only  needs  to  be  evaluated  once  for  each 
(■ m,n )  pair  and  then  stored  for  subsequent  use  in  evaluating  (61).  The  value  of  p  at  which  the  two 
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curves  are  forced  to  intersect  is  chosen  small  enough  that  the  small-p  approximation  is  accurate, 
p  =  0.1/n.  Evaluating  the  average  stopping  time  in  this  manner  forces  accuracy  as  p  — >  1  and  also 
for  p  <  0.1/n,  but  does  not  guarantee  accuracy  in  the  intervening  region.  For  example,  the  relative 
errors  for  the  n  =  7  case  shown  in  Figs.  3  and  4  have  a  maximum  error  of  e/j  =  0.18  when  p  =  0.4 
for  m  =  6.  While  an  18%  error  is  not  small,  the  error  is  often  less  and  there  is  little  alternative 
when  n  is  large  and  the  FSMP  exact  solution  is  not  easily  computed.  For  n  G  [6, 11],  the  errors  were 
generally  worst  for  m  =  n  —  1  (excepting  n  =  10  for  which  m  =  8  had  the  worst  error)  and  occurred 
for  some  p  G  (0.3, 0.5).  For  n  <  5,  the  worst  error  was  less  than  10%  and  occurred  for  smaller 
values  of  m  ( m  =  2  for  n  =  4  and  5  and  m  =  3  for  n  =  3).  Of  the  cases  evaluated  (n  G  [3, 11]),  the 
errors  increased  with  n. 

Of  the  techniques  available  for  approximating  the  average  stopping  time,  the  alternative  compu¬ 
tational  approximation  provides  the  best  overall  accuracy,  but  especially  so  for  values  of  p  >  0.1 
where  D  is  typically  evaluated.  For  evaluating  T  when  p  is  small,  the  small-p  approximation  should 
be  adequate;  however,  care  must  be  taken  if  n  is  large.  While  the  trend  toward  increasing  errors 
with  n  is  troubling,  there  is  little  alternative  for  evaluating  D  when  n  >  10. 


3.4  Choosing  an  optimal  m-of-n  detector 

The  standard  technique  for  detector  design  is  to  first  set  the  false-alarm  rate  (FAR)  and  then 
optimize  for  detection  performance.  However,  the  Neyman-Pearson  lemma,  which  dictates  that  a 
likelihood  ratio  yields  the  optimal  detector  under  these  conditions,  does  not  apply  to  the  sliding 
m-of-n  detector.  Optimality  will  still  be  defined  as  the  (m,  n)  pair  that  maximizes  detection 
performance  while  constraining  FAR.  Detection  performance  is  quantified  by  the  D  for  signals 
that  start  and  do  not  stop  or  have  very  long  durations  and  by  P, 1  for  finite-duration  signals. 


3.4.1  The  choice  of  m  dominates  the  FAR 


Typically,  p  will  be  small  under  the  null  hypothesis,  so  using  the  approximation  of  (56)  will  normally 
be  adequate.  Consider  the  logarithm  of  T  using  (56)  with  p  =  po, 

'  n  —  1 


log  T  ps  —m  logpo  —  (n  —  m  +  1)  log  (1  —  p o)  —  log 


m  —  1 


(65) 


If  po  is  small  enough,  the  first  term  dominates  the  second.  Then,  by  noting  that  the  last  term  is 
negative,  a  lower  bound  can  be  constructed  for  the  values  of  m  meeting  the  FAR  specification 


m  >  mo 


log  T 

logioT 

log  Po 

-logioPo 

(66) 
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where  [•]  is  the  ceiling  function.  For  a  fixed  value  of  m,  and  assuming  po  is  small  enough  that 
the  third  term  in  (65)  dominates  the  second,  increasing  n  decreases  T.  When  a  signal  is  present, 
increasing  n  (as  long  as  the  signal  duration  exceeds  n)  increases  the  opportunity  of  observing  m 
successes  and  therefore  increases  Pd-  As  previously  argued,  increasing  n  also  reduces  O.  Thus, 
the  optimal  value  of  n  for  a  fixed  rn  is  the  largest  value  still  meeting  the  false-alarm  specification. 
Define  this  value  as  n{m).  Noting  the  inequality  in  (66),  there  are  multiple  ( m,n )  pairs  that  will 
meet  the  FAR  specification  and  provide  the  largest  Pd  and  smallest  D  as  a  function  of  n  for  each 
associated  value  of  m.  The  following  two  sections  describe  how  to  choose  the  (m,  n)  pair  optimizing 
D  or  Pd. 


3.4.2  Minimizing  O 


Consider  the  example  FAR  specification  of  T  =  1020  with  po  =  10-3.  The  minimum  value  of  m 
satisfying  the  false  alarm  constraint  is  7,  easily  computed  from  |"log10(T)/ log10(po)]  •  As  seen  in 
Table  3,  the  largest  value  of  n  meeting  the  FAR  specification  for  m  =  7ish(7)  =  8.  Increasing  m 
results  in  successively  larger  values  of  the  optimal  n;  for  example,  h(10)  =  58  and  h(14)  =  292.  The 
larger  values  of  n  generally  result  in  less  overshoot  of  the  FAR  specification  and  therefore  might  be 
expected  to  provide  better  detection  performance;  however,  the  value  of  p\  appears  to  dominate  the 
choice,  as  seen  from  the  average  delay  to  detection  shown  in  Table  3.  When  p\  is  small,  (m,  n)  pairs 
with  larger  values  provide  the  best  performance  (e.g.,  when  p±  =  0.3,  D  is  minimized  by  m  =  10 
and  n  =  58).  However,  as  p\  approaches  one,  smaller  values  of  m  provide  the  best  performance. 
In  this  case,  the  optimality  of  the  (m,  n )  pair  with  the  smallest  value  of  m  satisfying  the  FAR 
specification  is  evident  by  noting  that  the  minimum  number  of  samples  to  a  decision  is  m. 

Table  3:  Pairs  of  (■ m,n )  meeting  the  false  alarm  specification  T  =  1020  for  po  =  10-3,  the  FAR 
achieved  as  computed  through  the  alternative  computational  approximation,  and  D  for  several 
values  of  p\.  The  bolded  entries  denote  the  (m,n)  pair  minimizing  D(pi). 


rn 

n{m) 

iogio  T 

D(0.30) 

D(0.50) 

D(0.70) 

0(0.95) 

7 

8 

20.128 

1713.7 

100.2 

20.6 

7.7 

8 

16 

20.150 

136.5 

24.0 

12.2 

8.4 

9 

32 

20.065 

45.5 

18.7 

12.9 

9.5 

10 

58 

20.022 

35.9 

20.0 

14.3 

10.5 

11 

96 

19.996 

37.0 

22.0 

15.7 

11.6 

12 

147 

19.995 

40.0 

24.0 

17.1 

12.6 

13 

212 

20.005 

43.3 

26.0 

18.6 

13.7 

14 

292 

20.013 

46.7 

28.0 

20.0 

14.7 

Values  of  m .  minimizing  D  as  a  function  of  T  and  pi,  denoted  as  rn* .  are  listed  in  Tables  4-6  for, 
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respectively,  po  =  10— 2 ,  10-3,  and  10~4.  The  associated  value  of  n  is  simply  n(m*),  the  largest 
value  of  n  meeting  the  FAR  specification  for  rn* .  As  seen  in  Tables  4-6,  m*  tends  to  be  mo  for 
smaller  values  of  T  or  larger  values  of  p\ .  For  example,  with  po  =  10~3,  setting  T  <  105  resulted 
in  choosing  mo  for  p\  >  0.1. 


Table  4:  Value  of  m  minimizing 

D  for 

Po  = 

10“2 

and 

various 

T  and  p\ 

togio  T 

mo 

m 

*  for 

Pi  = 

0.95 

0.9 

0.8 

0.7 

0.6 

0.5 

0.4 

0.3 

0.2 

0.1 

5 

3 

3 

3 

3 

3 

3 

3 

3 

4 

4 

5 

10 

5 

5 

5 

6 

7 

7 

7 

8 

9 

11 

15 

15 

8 

8 

9 

9 

10 

11 

11 

12 

14 

16 

23 

20 

10 

11 

12 

12 

13 

14 

15 

16 

18 

21 

31 

25 

13 

14 

15 

15 

16 

17 

18 

20 

22 

27 

39 

Table  5:  Value  of  m  minimizing 

D  for 

Po  = 

10“3 

and 

various 

T  and  p\ 

iogio  T 

mo 

m 

*  for 

Pi  = 

0.95 

0.9 

0.8 

0.7 

0.6 

0.5 

0.4 

0.3 

0.2 

0.1 

5 

2 

2 

2 

2 

2 

2 

2 

2 

2 

2 

2 

10 

4 

4 

4 

4 

4 

4 

4 

4 

5 

5 

6 

15 

5 

5 

5 

6 

6 

6 

7 

7 

7 

8 

10 

20 

7 

7 

7 

8 

8 

9 

9 

10 

10 

11 

13 

25 

9 

9 

10 

10 

10 

11 

11 

12 

13 

14 

16 

Table  6:  Value  of  m  minimizing 

D  for 

Po  = 

10-4 

and 

various 

T  and  p\ 

togio  T 

mo 

m 

*  for 

Pi  = 

0.95 

0.9 

0.8 

0.7 

0.6 

0.5 

0.4 

0.3 

0.2 

0.1 

5 

2 

2 

2 

2 

2 

2 

2 

2 

2 

2 

2 

10 

3 

3 

3 

3 

3 

3 

3 

3 

3 

3 

4 

15 

4 

4 

4 

4 

4 

5 

5 

5 

5 

5 

6 

20 

5 

5 

5 

6 

6 

6 

6 

6 

7 

7 

8 

25 

7 

7 

7 

7 

7 

8 

8 

8 

9 

9 

10 

While  the  optimal  value  of  m  may  be  taken  from  the  data  in  Tables  4-6  if  the  desired  values  of 
Po,  p\ ,  and  T  have  been  evaluated,  computation  of  D  for  each  value  of  m  >  mo  via  (50)  or  (61) 
(the  latter  when  n  >  10)  until  a  minimum  is  found  is  cumbersome  and  needlessly  time  consuming. 
An  alternative  approximate  method  may  be  found  by  exploiting  the  technique  proposed  in  the 
following  section  for  choosing  m  to  maximize  Pci  for  a  signal  of  a  specific  length.  By  assuming 
the  signal  length  is  equal  to  the  average  delay  to  detection  achieved  by  Page’s  test  operating  on 
Bernoulli  data  with  no  mismatch  in  the  choice  of  bias,  an  approximate  value  can  be  obtained  as 
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a  starting  point  to  find  to* .  The  average  delay  to  detection  for  Page’s  test  can  be  approximated 
using  the  efficacy  (rj)  from  (8)  and  the  FAR  specification, 


DP 


log  T 

V 


(67) 


The  value  of  to  approximately  maximizing  for  a  signal  of  length  Dp  for  the  example  of  Table  3  is 
shown  in  Fig.  6  as  a  red  circle  on  the  curves  of  D  against  to  with  the  optimal  value  of  m  denoted  by 
a  triangle.  For  p±  =  0.3,  this  method  yielded  m*  ~  9  where  the  optimum  value  was  10.  There  was 
a  similar  error  for  p\  =  0.5,  but  no  error  for  0.7  or  0.95.  The  error  counts  over  all  cases  considered 
in  Tables  4-6  are  reported  in  Table  7  where  there  appears  to  be  a  skew  toward  underestimating  to* 
by  one.  Noting  that  the  initialization  was  always  within  one  of  the  optimal  value,  there  are  clear 
reductions  in  computational  effort  compared  with  the  brute- force  search  when  m*  3>  mo-  As  seen 
in  Fig.  6,  it  is  still  important  to  find  the  optimal  value  of  m  definitively  as  erring  toward  a  smaller 
value  can  come  at  a  significant  penalty  in  D. 


Table  7:  Counts  of  occurrences  of  errors  in  the  initialization  approach  for  m*  over  all  the  cases 
shown  in  Tables  4-6. 


Po 

TO* 

-1 

0 

1 

1(T2 

2 

13 

35 

10~3 

1 

31 

18 

1(T4 

0 

43 

7 

Noting  the  rapid  rate  at  which  n(m)  increases  with  m  (e.g.,  see  Table  3),  only  the  technique  of  (61) 
is  amenable  to  extensive  evaluation  of  D.  However,  it  is  also  important  to  recall  that  (65)  is  an 
approximation,  requiring  that  po  be  small.  When  n(m )  is  large,  the  approximation  can  degrade  to 
the  point  where  it  is  not  accurate.  In  these  situations,  the  numerical  technique  described  by  (61) 
should  be  used  to  evaluate  T  as  well  as  D.  For  example,  the  column  in  Table  3  listing  the  FAR 
specification  is  computed  using  the  alternative  computational  approximation  of  (61)  for  the  ( m ,  n) 
pair  obtained  using  the  small-p  approximation.  While  the  inaccuracy  of  the  latter  is  small  in  this 
case,  it  is  still  evident  for  to  =11  and  12  and  can  be  worse  for  larger  values  of  po  or  n.  To  speed  up 
the  evaluation  of  n(m),  a  Newton- Raphson  iteration  can  be  applied  to  (65)  to  invert  the  equation 
for  n  by  writing  the  third  term  using  gamma  functions, 

logT  ~  — to.  log  po  —  (n  —  to  +  1)  log  (1  —  po)  —  logT  (to)  +  logT  (n  —  to  +  1)  —  logT  (n)  (68) 

and  noting  that  the  derivative  of  logr(.x)  is  the  digamma  function  ip(x). 
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Figure  6:  Average  delay  before  detection  as  a  function  of  m  for  po  =  10  3,  log10T  =  20,  and 
pi  =  0.3,  0.5,  0.7,  and  0.95.  The  optimum  value  of  m  increases  as  p\  decreases. 
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3.4.3  Approximately  maximizing  Pd 


Often  signals  persist  only  for  a  short  period  of  time,  in  which  case  the  probability  of  detecting  the 
signal  is  a  more  appropriate  performance  measure  than  the  average  delay  to  detection.  Consider 
the  example  of  the  previous  section  (recall  po  =  1CU3  and  log10  T  =  20)  with  an  individual-trial 
success  probability  of  p\  =  0.75.  As  noted  in  Table  3,  the  smallest  value  of  m  meeting  the  FAR 
specification  is  m o  =  7.  The  probability  of  detecting  a  signal  of  length  k  is  shown  in  Fig.  7  for  the 
7-of-7  and  7-of-8  detectors  (solid  lines).  The  dashed  curve  labeled  ‘7-of -k  envelope’  is  Pd  when  n 
equals  the  signal  length  and  represents  an  upper  bound  on  the  performance  of  all  detectors  with 
m  =  7  and  n  <  k.  Similarly,  the  curve  labeled  ‘8-of-A;  envelope’  is  an  upper  bound  when  m  =  8 
and  n  <  k.  While  Pd  for  rn-oi-n  detectors  where  k  >  n  must  be  evaluated  through  (40),  the  rn-oi-k 
performance  envelope  can  be  obtained  through  (43),  without  encountering  the  numerical  issues  of 
the  FSMP  approach  for  large  n.  This  is  exploited  in  developing  a  technique  for  choosing  the  (to,  n) 
pair  approximately  maximizing  Pd. 

In  Fig.  7,  it  can  be  seen  that  the  7-of-8  detector  provides  the  highest  Pd  for  signals  with  lengths 
k  <  9.  While  a  7-of-9  detector  would  have  higher  Pd  when  k  =  9,  h( 7)  =  8  (from  Table  3)  indicating 
that  the  7-of-9  detector  does  not  meet  the  FAR  specification.  When  k  =  10,  the  8-of-A;  envelope 
shows  that  the  8-of-10  detector  outperforms  the  7-of-8  detector.  Noting  that  h(8)  =  16,  the  8-of-A; 
detector  would  then  maximize  Pd  for  signals  of  length  10  <  k  <  16. 

Notice  that  when  the  signal  has  length  k  =  7,  both  the  7-of- 7  and  7-of-8  detectors  yield  the  same 
Pd  performance;  however,  they  have  different  FAR.  As  previously  noted,  increasing  n  when  to  is 
fixed  reduces  T.  As  such,  any  equality  in  Pd  for  detectors  with  the  same  value  of  to.  can  be  resolved 
by  choosing  the  smallest  value  of  n;  that  is,  the  smaller  of  the  signal  length  or  n(m). 

The  structure  of  the  Pd  curves  in  Fig.  7  provides  a  means  for  approximately  choosing  the  optimal 
(to,  n)  pair  to  maximize  Pd  while  constraining  T  without  resorting  to  evaluating  Pd  through  the 
FSMP  solution  found  in  (40).  Divide  the  signal- length  space  (the  natural  numbers)  into  consecutive, 
non-overlapping  intervals  C\,  £2,  ■  ■  ■ ,  such  that  within  any  interval  a  fixed  value  of  to  approximately 
maximizes  Pd-  Call  this  value  to*  for  the  zth  interval  and  let  the  associated  value  of  n  be 

n*  (to.,  k)  =  min  {A;,  h(m)}  .  (69) 

The  first  interval  is  C\  =  [1,  fcj]  with  to*  =  too  from  (66)  as  the  optimal  value  of  to.  The  next 
interval  starts  at  A;*  +  1  and  is  the  point  at  which  a  larger  value  of  to  has  a  performance  envelope 
exceeding  that  achieved  with  the  current  value.  In  Fig.  7,  the  first  interval  would  be  C\  =  [1,9] 
with  m\  =  7  and  n*( 7,  k).  The  second  starts  at  k\  + 1  =  10  with  =  8  and  n*(8,  k )  and,  although 
not  shown,  ends  at  k\  =  17.  The  example  is  continued  in  Fig.  8,  but  with  a  smaller  single-trial 
success  probability  p\  =  0.4  to  further  illustrate  the  optimization  process.  The  smaller  value  of  p\ 
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Figure  7:  Pri  as  a  function  of  signal  length  for  various  ( m,n )  pairs.  The  figure  illustrates  that  a 
7-of-n*(7,  k)  rule  works  best  when  the  signal  length  k  is  less  than  10  where  an  8-of-n*(8,  k)  rule 
then  works  best. 
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doesn’t  change  the  first  interval;  however,  the  second  interval  is  slightly  larger  at  £2  =  [10, 18]. 

Generalizing  the  description  of  the  signal-length  intervals,  the  ith  interval  ends  at  the  shortest 
signal  length  k*  such  that  a  larger  value  of  m  in  the  next  region  m*+1  >  m*  has  a  Pd  envelope 
greater  than  that  achieved  by  m*  at  k*  +  1  and  where  the  associated  value  of  n  meets  the  FAR 
specification  (i.e.,  k*  +  1  <  h  (m*+ , ) ) .  The  break-point  between  the  ith  and  (i  +  l)st  interval  can 
be  described  analytically  as 

k*  =  arg min  {k  :  Pd  (m*+1 ,  n*  (m -+1 ,  k))  >  Pd  (m* ,  n*  (rn* ,  k) ) }  -  1  (70) 

k 

where  the  Pd  envelope  is  (following  (43))  defined  as 

Pd(m,k)  =  Y,  („- lV(1  (71) 

j=m  ^  ' 


for  k  >  m  and  zero  otherwise. 

The  approximation  comes  from  using  the  Pd  envelope  (which  is  easily  computed)  held  constant  at 
h(m)  rather  than  the  exact  Pd  for  cases  where  k  >  n(rn).  In  the  example  of  Fig.  7,  the  interval  is 
correctly  chosen  resulting  in  an  optimal  selection;  however,  if  the  7-of-8  detector  performed  slightly 
better  at  k  =  10  than  the  8-of-/c  envelope,  the  rule  of  (70)  would  have  incorrectly  opted  for  the 
8-of-10  detector.  As  seen  in  Fig.  7,  the  Pd  envelope  has  a  higher  rate  of  ascent  than  the  exact 
Pd  for  cases  where  k  >  n(m),  so  the  sub-optimality  should  be  limited  to  the  interval  edges.  For 
applications  where  n  is  always  small  (e.g.,  n  <  10),  the  exact  Pd  can  be  computed  and  used  in 
place  of  the  Pd  envelope  on  the  right-hand-side  of  the  inequality  in  (70).  When  n  is  larger,  there  is 
no  clear  alternative  to  the  proposed  technique. 


UNCLASSIFIED 


approved  for  public  release 


31 


CausaSci  LLC 


UNCLASSIFIED 


2011-01 


Figure  8:  Pd  performance  envelopes  for  m-of-n  detection  as  a  function  of  signal  length  for  po  =  10~3, 
pi  =  0.4,  and  log10  T  =  20.  As  signal  length  increases,  the  optimal  value  of  m  increases  with  n  the 
lesser  of  the  signal  length  and  the  largest  value  of  n  meeting  the  FAR  specification  for  the  given  m. 
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4  Detection  performance  comparison 


The  detection  performance  of  m-of-n  detectors  is  compared  with  Page’s  test  both  in  terms  of  D 
and  Pd  for  a  given  FAR  specification.  It  has  been  proven  [6,  7]  that  Page’s  test  is  optimal  in 
terms  of  providing  the  lowest  D  for  a  given  T,  so  the  primary  issue  is  how  poorly  the  m-of-n 
detectors  perform  in  comparison.  As  there  is  no  such  optimality  with  respect  to  Pd,  an  analysis  of 
performance  is  required  to  determine  which  detector  performs  best  in  any  given  situation. 


4.1  Average  delay  to  detection 

ROC  curves  similar  to  those  found  in  Fig.  1  are  shown  in  Fig.  9  comparing  the  average  delay 
before  detection  as  a  function  of  the  average  time  between  false  alarms  for  the  optimal  (in  the 
sense  of  minimizing  D  for  a  constrained  T)  m-of-n  detector  and  Page’s  test.  The  case  of  po  =  10-3 
is  considered  for  various  values  of  p\ .  As  expected  based  on  its  optimality,  Page’s  test  is  seen  to 
provide  lower  D  than  the  best  m-of-n  detector  for  all  cases  considered,  with  an  increasing  disparity 
in  performance  as  p\  decreases. 

Both  detectors  require  an  assumption  on  signal  characteristics  (i.e.,  pi)  in  order  to  choose  the 
optimal  m-of-n  rule  or  the  bias  in  Page’s  test.  The  consequence  of  mismatch  as  a  function  of  the 
assumed  p\  is  shown  in  Fig.  10  for  po  =  10-3,  p\  =  0.5  and  various  FAR  specifications.  When 
the  design  value  of  p\  equals  the  true  value  (0.5  in  this  example),  the  performance  is  optimized; 
however,  mismatch  can  or  will  result  in  a  longer  average  delay  to  detection.  The  m-of-n  detectors 
resulted  in  zero  loss  in  many  of  the  cases  considered,  reflecting  the  quantization  involved  in  the 
design;  that  is,  the  same  values  of  m  and  n  will  be  optimal  for  a  range  of  conditions  rather  than 
a  single  point  condition.  Page’s  test  exhibited  an  immediate  performance  loss  upon  any  level  of 
mismatch,  but  still  outperformed  the  m-of-n  detectors  for  the  majority  of  cases  considered.  For 
both  detectors,  the  loss  arising  from  mismatch  appears  exacerbated  by  larger  values  of  T. 

From  the  perspective  of  the  average  delay  before  detection,  Page’s  test  is  clearly  a  more  desirable 
detector  than  the  sliding  m-of-n  detector.  It  enjoys  optimality  providing  the  lowest  D  of  all  tests 
meeting  the  FAR  specification,  but  also  exhibited  less  overall  sensitivity  to  mismatch  compared 
with  the  best  m-of-n  detector.  Additionally,  it  is  easier  to  implement,  only  requiring  a  simple 
recursion,  and  to  optimize.  Changing  the  FAR  specification  or  the  signal-presence  probability  only 
requires  a  change  in,  respectively,  the  detector  threshold  or  bias.  Changing  either  for  the  m-of-n 
detector  requires  a  completely  new  evaluation  of  the  optimum  values  of  m  and  n. 
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Figure  9:  ROC  curves  of  D  as  a  function  of  log10  T  for  the  best  m-of-n  detector  and  Page’s  test 
for  po  =  10“3  and  various  values  of  p\ . 
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Figure  10:  D  under  a  mismatch  between  assumed  and  actual  values  of  p\  for  the  best  m-of-n 
detector  and  Page’s  test  for  po  =  10~3,  various  values  of  log10  T  and  a  true  value  of  p\  =  0.5. 
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4.2  Probability  of  detection 


Based  on  the  optimality  of  Page’s  test  for  minimizing  D  and  the  Neyman-Pearson  optimality  of  the 
likelihood  ratio  test  for  detection  with  fixed  sample  sizes,  it  is  easy  to  assume  that  the  likelihood 
ratio  would  be  at  the  heart  of  the  detector  optimizing  Pd  for  a  finite-duration  signal  in  a  sequential 
detection  setting  where  signal  onset  is  unknown.  Interestingly,  the  m-oi-n  detectors  appear  to 
provide  comparable  or  better  performance  than  Page’s  test. 


As  described  in  Sect.  3.4.3,  evaluating  Pd  for  m-oi-n  detectors  is  difficult  when  n  is  above  about 
10,  which  forces  the  use  of  the  Pd  envelope  described  by  (70)  with  the  breakpoints  in  (71)  dictating 
when  to  change  m  and  n  chosen  from  (69).  The  resulting  function  is  a  lower  bound  on  Pd  as  a 
function  of  signal  length.  In  a  practical  sense,  Pd  above  this  lower  bound  can  only  be  obtained 
when  the  signal  length  is  small  enough  to  allow  exact  evaluation  of  Pd  through  the  FSMP  approach 
or  if  simulation  analysis  is  undertaken  to  approximate  Pd- 


For  Page’s  test,  Pd  is  evaluated  using  the  FSMP  approach  described  in  [8]  applied  to  the  binary 
data  and  data  mapping  function  of  (18).  The  approach  requires  quantizing  the  Page’s  test  detection 
statistic  (Y);  in  (2))  and  data  mapping  function  in  order  to  characterize  Page’s  test  as  an  FSMP, 
similar  to  what  was  done  for  the  m-oi-n  detector  in  Sect.  3.1.  The  data  mapping  function  of  (18) 
can  only  produce  values  of  1  —  b(pQ.p\ )  and  —  b(po,pi),  so  the  transition  matrix  for  the  continuing 
states  is  easily  constructed.  However,  accuracy  in  the  result  requires  a  fine  enough  quantization  for 
there  to  be  many  (at  least  10)  steps  between  either  value  and  zero.  This  limits  the  practical  utility 
of  the  method  for  cases  where  b(po,pi)  is  near  one  or  zero.  In  evaluating  Pd  for  both  detectors, 
so  called  ‘latent’  detections  (those  occurring  after  the  signal  has  stopped,  but  while  the  detection 
statistic  is  still  influenced  by  signal  presence)  are  ignored. 


Consider  the  case  of  po  =  10~3  and  log10  T  =  10  shown  in  Fig.  11  where  Pd  or  its  lower  bound 
are  shown  as  a  function  of  signal  length  for  several  values  of  p\.  Generally,  the  m-oi-n  detectors 
provide  equivalent  or  better  detection  performance  than  Page’s  test.  Pd  for  Page’s  test  appears  to 
initially  track  the  Pd  envelope  for  the  m-of-n  detectors  where  m  is  chosen  as  the  smallest  number 
of  Page-test  updates  allowing  a  detection  declaration, 


mp  = 


h 

1  -  b(p0,pi) 


(72) 


However,  it  shifts  to  follow  the  Pd  envelope  for  the  next  higher  value  of  m  before  the  breakpoint 
technique  described  in  Sect.  3.4.3  for  the  m-oi-n  detectors.  As  seen  in  Fig.  11,  the  m-oi-n  detector 
can  provide  a  significant  performance  improvement  over  Page’s  test  for  signals  with  lengths  in 
these  regions.  The  Pd  curves  for  pi  =  0.5  are  expanded  in  Fig.  12  and  shown  along  with  the  4-of-A: 
detector  curves  up  to  k  =  10.  This  figure  illustrates  that  the  Pd  envelope  is  a  lower  bound  for 
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performance,  with  the  4-of-10  detector  providing  better  performance  than  the  5-of -k  (which  would 
be  chosen  by  the  Pt 1  envelope  approach)  performance  when  k  =12  through  15.  Thus,  the  range  of 
signal  lengths  over  which  the  m-of-n  detectors  outperform  Page’s  test  is  greater  than  that  indicated 
by  the  Pd  envelope,  but  not  necessarily  attainable  if  it  requires  evaluating  Pd  through  the  FSMP 
approach. 

It  can  be  seen  in  Fig.  12  that  if  the  signal  has  length  6,  Page’s  test  and  the  4-of -k  detectors  for 
k  =  4, . . . ,  10  all  provide  the  same  Pd-  As  noted  in  Sect.  3.4.3,  choosing  the  4-of-4  detector  results 
in  the  best  achieved  FAR.  This  can  be  seen  in  Fig.  13  where  it  clearly  has  the  highest  average  time 
between  false  alarms  not  only  for  the  design  value  of  po  =  1CU3,  but  also  if  the  observed  value  of 
Po  differs  from  the  design  value.  In  fact,  the  4-of-4  detector  is  robust  enough  to  meet  the  FAR 
specification  even  when  the  observed  po  is  three  times  greater  than  the  design  value.  As  seen  in 
Fig.  13,  the  robustness  decreases  as  n  — >  n(m),  the  maximum  value  meeting  the  FAR  specification 
for  the  given  value  of  to.  However,  from  Figs.  11  and  12,  the  m-of-n  detectors  have  higher  Pd  than 
Page’s  test  in  this  region  of  signal  length. 

The  primary  disadvantage  of  the  m-of-n  detectors  with  respect  to  maximizing  Pd  is  that  different 
(m,  n)  pairs  are  required  to  optimize  the  performance  for  each  signal  length.  The  performance 
results  for  Page’s  test  seen  in  Figs.  11-13  are  obtained  with  just  one  design  (i.e.,  bias  and  threshold) 
for  each  value  of  pi,  irrespective  of  the  signal  length.  Noting  the  significant  loss  in  performance  seen 
in  Fig.  12  when  n  is  less  than  the  signal  length,  it  is  better  to  overestimate  the  signal  length  if  it  is 
not  known  precisely  or  to  use  Page’s  test.  However,  it  is  important  to  maintain  the  optimal  value 
of  to;  that  is,  overestimating  to  the  point  where  to  increases  will  result  in  a  lower  Pd  while  merely 
increasing  n  such  that  n  <  n(m)  will  maintain  Pd  and  trade  robustness  in  FAR  for  Pd  optimality 
over  a  broader  range  of  signals  (e.g.,  using  the  4-of-10  detector  optimizes  Pd  for  signals  of  length 
k  <  10  and  just  meets  the  FAR  specification).  Thus,  there  are  clear  advantages  in  using  the  m-of-n 
detectors  rather  than  Page’s  test  when  the  signal  length  is  known. 
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Figure  11:  Pd  as  a  function  of  signal  length  for  Page’s  test  and  the  Prj  envelope  (lower  bound)  for 
the  optimal  m-of-n  detectors  for  po  =  10— 3 ,  log10T  =  10  and  various  values  of  p±. 
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Figure  12:  Pd  as  a  function  of  signal  length  for  Page’s  test  and  various  m-of-n  detectors  for  po  = 
1CT3,  log10  T  =  10  and  pi  =  0.5. 
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Figure  13:  FAR  achieved  with  detectors  designed  for  po  =  10~3,  log10T  =  10  and  p\  =  0.5  when 
the  observed  value  of  po  differs  by  up  to  an  order  of  magnitude.  Some  of  the  m-of-n  detectors  are 
robust  enough  to  meet  the  FAR  specification  even  when  the  observed  po  is  over  three  times  larger 
than  the  design  value. 
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5  Conclusions 


This  report  has  presented  analysis  techniques  for  the  sliding  m-of-n  detector  to  enable  design  meet¬ 
ing  a  FAR  specification  and  approximately  optimizing  detection  performance.  Exact  computational 
techniques  were  presented  for  evaluation  of  D,  T,  and  Pd  through  the  FSMP  description  of  the 
sliding  m.-of-n  detector  presented  in  [3]  under  stationary  data  conditions.  As  the  number  of  states 
grows  exponentially  with  n,  these  techniques  become  computationally  cumbersome  for  n  >  10  and 
impractical  shortly  thereafter.  Approximations  from  scan-statistics  literature  [5]  may  be  useful  for 
larger  values  of  n,  but  are  still  limited  at  smaller  values  of  p  and  assume  a  minimum  detection 
time  of  n  rather  than  m,  restricting  their  utility  to  evaluation  of  T  for  moderate  values  of  p,  un¬ 
less  extensions  exist  to  accommodate  this  alternative  stopping  criteria.  To  counter  the  numerical 
limitations  of  these  approaches,  approximations  applicable  for  a  much  larger  range  of  n  were  devel¬ 
oped,  including  a  small-p  approximation  useful  for  evaluating  T  and  an  alternative  computational 
approximation  for  D. 

The  value  of  m  was  seen  to  dominate  the  FAR  of  the  sliding  m-of-n  detector,  although  multiple 
(m,  n)  pairs  with  m  above  a  minimum  value  were  found  to  meet  the  specification.  Of  these,  a 
unique  pair  can  be  chosen  to  either  minimize  D  or  approximately  maximize  Pd-  The  latter  was 
done  by  choosing  the  (m,  n )  pair  maximizing  an  easily  evaluated  lower  bound  on  the  probability  of 
detecting  a  finite-duration  signal  using  m.-of-n  detectors.  Given  the  optimal  m,  any  value  of  n  from 
the  signal  length  up  to  h(m)  (the  maximum  value  meeting  the  FAR  specification)  maximizes  the 
Pd  lower  bound,  with  the  lower  values  having  greater  FAR  robustness  but  a  smaller  range  of  signal 
lengths  over  which  they  are  optimal.  Designing  the  detectors  to  minimize  D  was  less  complicated 
with  the  optimal  value  of  m  approximated  by  that  optimizing  Pd  for  a  signal  of  length  equal  to 
the  average  delay  to  detection  achieved  by  Page’s  test.  The  value  of  n  minimizing  D  was  then  the 
largest  one  meeting  the  FAR  specification. 

While  the  m-of-n  detectors  operate  on  binary  or  decision-level  data,  Page’s  test  can  be  applied  to 
binary  data  or  the  pre-thresholded  data.  An  example  illustrated  the  significant  loss  in  performance 
when  decision-level  data  are  used  unnecessarily.  The  performance  comparison  between  Page’s  test 
operating  on  binary  data  and  the  sliding  m-of-n  detectors  showed  that  the  latter  can  come  close 
to  the  optimal  performance  of  Page’s  test  with  respect  to  the  average  delay  to  detection,  but  not 
exceed  it  except  under  extreme  mismatch.  Although  the  sliding  m-of-n  detector  was  robust  to  minor 
mismatch,  often  yielding  no  loss,  it  could  result  in  significant  loss  at  larger  levels  of  mismatch.  The 
sliding  m-of-n  detectors  clearly  outperformed  Page’s  test  with  respect  to  the  probability  of  detecting 
finite-duration  signals.  In  cases  where  the  detection  performance  was  similar,  the  sliding  m-of-n 
detectors  exhibited  significant  FAR  robustness  to  varying  individual  trial  false-alarm  probability 
(po)  where  Page’s  test  did  not.  Although  not  considered  extensively  in  this  report,  design  of  the 
threshold  used  to  produce  decision-level  data  should  be  incorporated  into  the  optimization  process. 
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An  initial  approximation  was  obtained  by  approximately  optimizing  the  asymptotic  efficiency  of 
Page’s  test. 

The  primary  disadvantage  of  the  sliding  m-of-n  detectors  is  that  optimal  design  requires  different 
configurations  (i.e.,  (m,  n)  pairs)  that  depend  on  the  FAR  specification,  the  individual-trial  suc¬ 
cess  probabilities  and  signal  length.  With  Page’s  test’s  optimality  for  minimizing  D  and  simpler 
design,  it  remains  the  most  desirable  detector  for  signals  of  unknown  or  infinite  length.  However, 
the  performance  improvement  attained  by  a  properly  designed  m-of-n  detector  over  Page’s  test 
for  detecting  finite-duration  signals  was  significant  enough  that  they  are  the  clear  choice  when 
signal  length  is  known  or  known  within  a  range  where  one  optimal  sliding  m-of-n  detector  can  be 
implemented. 
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A  MGF  unity  roots 

A-l  Exponential  data 


Finding  the  unity  root  of  the  MGF  is  equivalent  to  finding  the  root  (zero)  of  the  logarithm  of  the 
MGF,  which  is  the  cumulant  generating  function  (CGF).  For  the  exponentially  distributed  data 
described  in  Sect.  2,  the  CGF  of  g(X )  under  the  signal-present  hypothesis  is 

Ts(t)  =  -b(s)t  -  log  [1  -  A  (1  +  s)t\  (73) 

for  t  <  [A(l  +  s)]"1.  The  CGF  under  the  null  hypothesis  is  simply  To (t),  which  still  depends  on 
the  design  SNR  s.  Under  the  null  hypothesis  and  the  condition  stated  in  (4),  Eq  [(7(A)]  <  0,  it  can 
be  shown  that  to  >  0.  Similarly,  1 1  <  0  under  the  alternative  hypothesis  when  E\  [5(A)]  >  0.  Note 
that  some  values  of  s  <  s  may  violate  the  condition  stated  in  (4)  yielding  t\  <  0,  and  very  poor 
detection  performance. 


The  procedure  for  obtaining  t\  is  described  below;  to  may  be  obtained  by  setting  s  =  0  without 
altering  s.  The  root  ti  (Ts(ti)  =  0)  is  easily  calculated  through  a  Newton-Raphson  iteration  (NRI), 


1 1  1\ 


Ts(ti) 


(74) 


where  the  the  derivative  of  the  CGF  is 

<(t)  =  ~b(S)  + 


[A  (l  +  ^-'-t' 


(75) 


The  NRI  must  be  initialized  far  enough  away  from  the  origin  that  it  converges  on  the  non-zero  root 
as  opposed  to  the  root  at  the  origin.  This  is  accomplished  by  starting  beyond  the  minimum  point 
of  the  CGF  occurring  at 


1 


1 


(*  =  a(TT7)  "  h5> 

which  sets  the  derivative  in  (75)  to  zero.  When  the  root  is  negative,  the  NRI  can  be  started  at  2f*; 
however,  when  it  is  positive,  it  must  be  started  between  t*  and  the  upper  limit  on  the  range  of  t 
for  which  the  MGF  is  defined. 


A-2  Bernoulli  data 


The  CGF  of  the  Bernoulli  data  with  probability  p  is 

Tp(f)  =  log  |pet(1~^  +  (1  -  p)e~tb 
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and  its  derivative  is 


K(t)  = 


pjl  -  -  (1  -p)be~tl 

pei(l-b)  _|_  ^  _  pje-tb 


(78) 


The  non-zero  roots  are  again  found  using  a  NRI  as  described  in  (74).  When  the  root  tp  >  0  (i.e., 
p  <  b),  the  iteration  can  be  initialized  using  the  approximation 


pe 


^p(tp)  ~  log 
Equating  (79)  to  zero  and  solving  results  in 

-  log  (p) 


tp(l-b) 


tp  — 


1-6 


Alternatively,  when  tp  <  0  (i.e.,  p  >  b)  the  approximation 


yields  the  NRI  initialization 


'l'p(tp)  «  log  (1  —  p)e  tpb 
log(l  -p) 


tp  — 


(79) 

(80) 

(81) 

(82) 
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B  Matlab  code  for  ra-of-n  detector  performance  evaluation  and 
design 


B-l  Subroutines  to  represent  the  sliding  rn-of-n  detector  as  an  FSMP 


ramnnrammra  nmm 


function  S=mofn_setup_opt (m,n) 

•/.  S=mofn_setup_opt (m,n) 

%  Returns  a  structure  S  of  states  characterizing  a  sliding  m-of-n  detector. 
%  S.ns  =  number  of  states 

°i  S.H  =  integer  representation  of  the  state 
*/.  S.Hp  =  transitioning  state  upon  a  success 
°/0  S.iHp  =  index  of  state  in  S.Hp 
°/0  S.Hm  =  transitioning  state  upon  a  failure 
1  S.iHm  =  index  of  state  in  S.Hm 

“/«  S.P  =  temporary  storage  for  Williams  (DST0-TN-0132)  algorithm  for  Pd 


%  Data  are  stored  in  a  .mat  file  so  only  have  to  do  this  once 
fname= [,mofn_struct_n)  int2str(n)] ; 
if  (exist ( [fname  ’ .mat1] , ’file’)==2) , 
load(fname) 
else 

%  Set  structures  to  be  empty  for  each  m  on  first  use 
for  i=l:n, 

Sm{i}=  []  ; 
end; 
end; 

'/.  Calculate  structure  if  not  already  stored 
if  isempty(Sm-(m})==l , 

dispC -  Stand  by,  calculating  structures  - 

fprintfC  m=,/0i  n=°/0i  1  ,m,n)  ; 

%  Add  the  accepting  state 

i=l; 

S ,ns=i ; 

S.P(i)=0; 

S ,H(i)=2~n; 

S .Hp(i)=2~n; 

S .Hm(i)=2~n; 

%  Number  of  non-accepting  states  (eq  2  in  DSTD-TN-0132) 

%  Easier  to  enumerate  and  find  number  of  non-accepting  states 
hiall=0 :  (2~  (n-1) -1)  ;  "/«  integer  rep  of  each  possible  history 
ibin=binnum(hiall)  ;  °/n  sum  of  l’s  of  binary  rep  of  each  state 
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hna=hiall (ibin<m) ;  7.  integer  equivalents  of  non-accepting  states 

*/.  Loop  through  all  possible  non-accepting-states 

nnas=length(hna) ; 

fprintfC#  states  =  7oi  \n',nnas); 

for  i=l:nnas, 

if  (nnas>500)M(rem(i,round(nnas/10))==l)  , 
fprintf 
end; 

7,  Reduce  to  smallest  equivalent  state 
hnaO=reduce_hexad(hna(i) ,m,n) ; 

/  If  it  hasn't  been  added,  then  add  it 
if  ~any (S ,H==hnaO) , 

S=addhexad_opt(S,hnaO,m,n) ; 
end; 
end; 

7.  Sort  hexads  by  H 
[~ , isrt] =sort (S .H) ; 

S .H=S .H(isrt) ; 

S .Hp=S .Hp(isrt) ; 

S .Hm=S .Hm(isrt) ; 

S.P=S.P(isrt) ; 

7o  Add  indices  for  each 
S . iH=l : S .ns ; 
for  i=l:S.ns, 

S.iHp(i)=find(S.Hp(i)==S.H) ; 

S. iHm(i)=find(S.Hm(i)==S.H) ; 
end; 

7»  Save  data  to  .mat  file  for  later  retrieval 
Sm-[m}=S ; 

eval(['save  ’  fname  ’  Sm']); 

7.  If  structure  is  stored,  then  return  it 
else 

S=Sm{m> ; 
end; 

[S . H 1  S.Hm’  S . Hp '  S.iHm’  S.iHp’]; 


function  S=addhexad_opt(S,h,m,n) 

7,  S=addhexad_opt  (S  ,h,m,n) 

7.  Function  internal  to  mofn_setup_opt (m,n)  to  add  a  hexad  to  structure 

7. 

i=S.ns+l;  7o  Find  index  for  this  hexad 
S.ns=i;  °L  Increment  number  of  hexads 
S.P(i)=0; 

S.H(i)=h;  7.  Current  hexad’s  integer  state  represenation 
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7  If  next  trial  is  a  failure 
hm=mod(2*h, 2~ (n-1) ) ;  7,  Integer  state 

hmO=reduce_hexad(hm,m,n)  ;  %  Reduce  to  smallest  equivalent  state 
S . Hm ( i ) =hmO ; 

7„  If  next  trial  is  a  success 
•/.  If  meeting  m-of-n  criteria 
if  binnum(h)==(m-l) , 

S.Hp(i)=2~n;  %  Accepting  state  integer  ID 

•/.  Otherwise  into  a  continuing  state 
else 

hp=hm+l;  %  Integer  state 

S .Hp(i)=reduce_hexad(hp,m,n) ;  7.  Reduce  to  smallest  equivalent  state 
end; 


7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.77.7.7.7.7.7.7.7.7.7.7. 


function  hO=reduce_hexad(h,m,n) 

7.  hO=reduce_hexad(h,m,n) 

7.  Function  internal  to  mofn_setup_opt (m,n)  to  reduce  a  state  to  the 
7.  "lowest  common  state",  which  minimizes  the  number  of  unique  states 
7.  in  the  m-of-n  FSMP 


hO=h; 

if  binnum(h) ~=(m-l) , 
hp=mod(2*h,2~ (n-1) )+l ; 
np=reduce_hexad(hp,m,n) ; 
if  np<h, 

hO=f loor (np/2) ; 
end; 
end; 


7.7.7.7.7.7.7.7.7.77.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7. 


function  b=binnum(i) 

7  b=binnum(i) 

7.  Returns  the  number  of  1  ’  s  in  the  binary 
7.  represetnation  of  integer  i 
7.  binnum(2)=l 
7.  binnum(3)=2 
7.  binnum(4)=2 

7.  Input  i  can  be  a  vector  (output  is  a  column  vector) 
x=dec2bin(i ( : ) ) ; 

[nr , nc] =size  (x)  ; 

b=sum(reshape(str2num(reshape(x,nr*nc,l)) ,nr,nc) ’) ’ ; 
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B-2  Subroutines  for  performance  measures  from  FSMP 


function  [uK , sK, err] =mof n_stats (pv, S , POv) 

•/.  [uK,sK]=mofn_stats(pv,S,POv) 

7„  Returns  the  mean  and  standard  deviation  of  the  stopping  time 
7„  for  the  m-of-n  process  defined  by  the  structure  S. 

%  pv  =  vector  of  Bernoulli  probabilities  to  be  evaluated 
7«  S  =  structure  describing  m-of-n  process  (see  mofn_setup) 

°/0  POv  =  vector  of  initial  starting  probabilities 
1  =  [1  0  ...  0] ’  by  default 

7. 


*/.  uK  =  average  stopping  time  (latency)  for  each  value  in  pv 
°/0  sK  =  standard  deviation  of  stopping  time 

l 


[Pep, Pern, pep] =mofn_states (S) ; 
ns=length(pcp)+l ; 
if  nargin<3, 

P0v=[l;zeros(ns-2,1)]  ; 
else 

P0v=P0v(l : (ns-1) ) ; 
end; 

uK=zeros(size(pv)) ; 
sK=uK ; 

for  i=l : length (pv) , 

7,  Construct  continuing  state  trans  matrix 
Pcc=pv(i)*Pcp+(l-pv(i) )*Pcm; 
pcs=pv(i)*pcp; 

7,  Eigen  decomposition 
[Vc,Dc]=eig(Pcc) ; 
dc=diag(Dc) ; 

70  If  eigenvectors  are  not  linearly  independent  use  pseudoinverse 
if  rcond(Vc)<le-10, 
cc=pinv (Vc) *pcs ; 
else 

cc=Vc\pcs ; 
end; 

7,  Check  to  make  sure  pseudoinverse  worked 
errl=max(abs(pcs-Vc*cc)) ; 

if  errl>sqrt (eps) ,  disp(’****  Error  ****’);  end; 
err (i)=rcond(Vc) ; 

%  Form  coefficients 
c=(P0v’*Vc) . ’ . *(cc) ; 

%  Form  mean,  power  &  standard  deviation 
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uK(i)=real (sum(c . / (1-dc) . ~2) ) ; 
pK=real(sum(c . * (1+dc) . / (1-dc) . ~3) )  ; 
sK(i)=sqrt (pK-uK(i) ~2 ) ; 
end; 


nniTaninninniniinni 


function  Pd=mofn_pd(pv,Lv,m,n) 

7  Pd=mofn_pd(pv,Lv,m,n) 

7„  Returns  Pd  for  set  of  success  probabilities  in  pv  for  the  time 
•/.  indexes  in  the  vector  Lv 

7  Output  Pd  is  a  matrix  with  length (pv)  rows  and  length (Lv)  columns 

7. 


S=mofn_setup_opt (m,n) ; 

7  Don’t  do  this  if  there  are  too  many  states 
if  S.ns>1010, 

Pd=nan*ones(size(Lv) ) ; 
else 

[Pep, Pern, pcp]=mofn_states(S) ; 
ns=length(pcp)+l ; 

P0v= [1 ; zeros (ns-2 , 1)] ; 
nL=length(Lv) ; 
np=length(pv) ; 

Pd=zeros(np,nL) ; 

7.  Construct  continuing  state  trans  matrix 
for  i=l:np, 

Pcc=pv(i) *Pcp+ (l-pv(i) ) *Pcm; 
pcs=pv(i)*pcp; 

7.  Eigen  decomposition 
[Vc ,Dc] =eig(Pcc) ; 
dc=diag(Dc) ; 

c=conj (Vc ’ *P0v) . * (Vc\pcs) ; 

Pd(i , : )=abs (sum(repmat (c . / (1-dc) , 1 ,nL) . * (1-repmat (dc , 1 ,nL) . ~repmat (Lv(:)’,ns-l,l)))); 
end; 
end; 


7.7.7.7.7.7.77.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7. 


function  [Pep , Pcm , pep] =mof n_states (S) 

7  [Pcp,Pcm,pcp]=mofn_states(S) 

7  Returns  the  indicator  matrices  and  vector  required  to  form  the 

7  FSMP  transition  probability  matrix. 

ns=S.ns; 

Pcp=zeros (ns-1) ; 

Pcm=Pcp; 

pcp=zeros (ns-1 , 1) ; 
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for  i=l : (ns-1) , 
if  S . iHp(i)==ns , 
pcp(i)=l ; 
else 

Pcp(i , S . iHp(i) )=1 ; 
end; 

Pcm(i , S . iHm(i) )=1 ; 
end; 


B-3  Subroutines  for  performance  measure  approximations 

function  [uK , sK] =mof n_st at s_asy (p , m , n) 
l  [uK, sK] =mofn_stats_asy(p,m,n) 

%  Returns  the  mean  and  standard  deviation  of  the  stopping  time 
%  for  an  m-of-n  process  when  the  probability  of  success  is  small 
%  (i.e.,  pv  ->  zero) 

1 

%  Primarily  useful  in  evaluating  average  time  between  false  alarms. 

1 

l  p  =  Bernoulli  probability  of  success 
V.  m  =  number  of  successes  required  for  stopping 
l  n  =  window  over  which  successes  are  counted 
°/«  Any  input  can  be  a  vector 
1 

%  uK  =  average  stopping  time  (latency)  for  each  value  in  pv 
%  sK  =  standard  deviation  of  stopping  time 
1 


%  Evaluate  (n-1) -choose- (m-1)  using  gamma  functions  to  allow  vector  inputs 
if  m==l, 

pa=p; 

else 

pa=exp(gammaln(n) -gammaln(m)-gammaln(n-m+l) ) . * (p . ~m) . * ( (1-p) . ~ (n-m+1) ) ; 
end; 


uK=l . . /pa; 
sK=sqrt (1-pa) ./pa; 


mmmmrarafflmmm 


function  [uK, sK] =mofn_stats_cf cn(p,m,n,delmin) 

'/.  [uK,  sK]  =mofn_stats_cf  cn(p,m,n) 

%  Returns  the  mean  and  standard  deviation  of  the  stopping  time 
%  for  an  m-of-n  process  when  the  probability  of  success  is  small 
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°/0  (i.e.,  pv  ->  zero) 

7„  p  =  Bernoulli  probability  of  success 
1  m  =  number  of  successes  required  for  stopping 
1  n  =  window  over  which  successes  are  counted 
%  delmin  =  minimum  step  size  to  evaluate  integral 
1 

%  uK  =  average  stopping  time  (latency)  for  each  value  in  p 

%  sK  =  standard  deviation  of  stopping  time 

°/. 

if  nargin<4, 
delmin=  []  ; 
end; 

if  m==l, 

pa=exp( (gammaln(n) -gammaln(m) -gammaln(n-m+l) ) ) . * ( (p . ~m) ) ; 
else 

pa=(l/m)*p.*exp(-(m-l)*mofn_cfcn_int(-log(p) ,mofn_cf cn(m,n) , delmin)) ; 
end; 

uK=l . . /pa; 
sK=sqrt (1-pa) ./pa; 


mmmmmmmmmm 


function  c=mofn_cf cn(mv,nv) 

•/.  c=mofn_cf cn(mv,nv) 

•/.  Returns  the  value  of  c(m,n)  for  each  element  in  the  vectors  mv  &  nv 


1  Make  vectors  the  same  size  if  one  is  a  scalar 
nmax=max(nv) ; 

if  (length (nv)==l)&& (length (mv) >1) 
nv=nv*ones (size (mv) ) ; 
elseif  (length(mv)==l)&&(length(nv)>l) , 
mv=mv*ones (size (nv) ) ; 
end; 

c=zeros (size (nv) ) ; 

•/.  Load  file  see  if  data  for  this  case  is  present 
fname=’mofn_cf  cn_data’  ; 
if  (exist ( [fname  1 .mat1] , ’file’)~=2) , 
for  i=l:nmax, 

Cs{i}=[0;  nan*ones(i-l , 1)] ; 
end; 
else 

load(fname) ; 
end; 

%  Populate  additional  cells  w/  nans  up  to  n 
nCs=length(Cs) ; 
for  i=(nCs+l) :nmax, 
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Cs-[i}=[0;  nan*ones (i-1 , 1)]  ; 
end; 

opts=optimset ( ’TolX’ , le-8) ; 

Iogl0e=logl0(exp(l)) ; 
delmin=0 . 01 ; 
qnew=0 ; 

for  i=l : length (nv) , 
n=nv ( i ) ; 
m=mv ( i ) ; 

7,  See  if  data  for  present  case  exists 
if  ~isnan(Cs{n}(m) ) , 
c(i)=Cs{n}(m) ; 
else 

qnew=l ; 

/  Get  best  fit  exponent 
p0v=0 . 1/n; 
u0v=log(p0v) ; 

7.  Use  asymptotic  approximation  as  barometer 

Ft=-logl0e* (gammaln(n) -gammaln(m) -gammaln(n-m+l) )-m*logl0(p0v) -(n-m+1) *logl0(l-p0v) ; 
/  Find  optimal  exponent 

c(i)=fminbnd(@(c) ((Iogl0(m)-logl0(p0v)+. . . 

(m-1) *mofn_cf cn_int (-uOv, c .delmin) *logl0e-Ft) . ~2) , 0 ,n-l , opts) ; 

Cs{n>(m)=c(i) ; 
end; 
end; 

7.  Save  the  data  if  any  new  ones  added 
if  qnew, 

eval([’save  ’  fname  ’  Cs’]); 
end; 


function  F=mofn_cf cn_int (xu, c ,delmin) 

7.  F=mofn_cf  cn_int  (xu,  c ,  delmin) 

7.  Subroutine  used  by  c=mofn_cf cn(mv,nv) 


7.  Spacing  of  approximations  in  -log(p) 
if  nargin<3, 
delmin=  []  ; 
end; 

if  "any (delmin) ,  delmin=0.25;  end; 
del=min( [delmin  c/2]); 

7.  Center  points  for  quadratic  Taylor  approximation 
xc=(del/2) :del: (max(xu)+del/2) ; 

[f c ,fplc , fp2c] =mofn_cf cn_deriv(xc , c) ; 

7.  Coefficients  of  quadratic  Taylor  approximation 
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b0=(fp2c/6)  ; 

bl=0 . 5* (fplc-xc . *fp2c) ; 

b2=f c-xc . *fplc+fp2c . * (xc . ~2) /2; 

7  Evaluation  points  for  integral  at  region  edges 
xim=xc-del/2 ; 
xip=xc+del/2 ; 

7„  Integral  over  each  region 

Fpm=bO . * (xip . ~3-xim. ~3)+bl . * (xip . “2-xim. ~2)+b2 . * (xip-xim) ; 

7  Find  nearest  point  for  each  element  of  xu 
F=zeros (size (xu) ) ; 
for  i=l : length (xu) , 

[~ , iu] =min(abs (xip-xu(i) ) ) ; 

7.  Integral  from  last  region  edge  to  desired  point 

Fu=bO (iu) * (xu(i) ~3-xip(iu) ~3)+bl (iu) . * (xu(i) ~2-xip(iu) ~2)+b2(iu) . * (xu(i) -xip(iu) ) ; 
7o  Final  result 
F(i)=sum(Fpm(l : iu) )+Fu; 
end; 


mnny.nmnmy.mnmm 


function  [f , f pi , f p2] =mof n_cf cn_deriv (x , c) 
7.  [f  ,  fpl ,  fp2]  =mofn_cf  cn_deriv(x,  c) 

7  Function  used  by  mofn_cf cn_int 

ex=exp(x) ; 

f=(l-exp(-x) ) . ~c ; 

fpl=c*f . / (ex-1) ; 

fp2=fpl . * (c-ex) . / (ex-1) ; 


7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.77.7.7.7. 


function  [pd , k] =mof n_pd_k_less_n (p , m , n) 

7,  [pd,k]  =mofn_pd_k_less_n(p,m,n) 

7  Returns  Pd  for  the  case  of  k<=n  using  exact  formula 
7  p  =  scalar  probability  of  success 

7  k  =  1  :n 


k=l :n; 

pd=zeros(size(k) ) ; 
kl=m:n; 

pd(kl)=cumsum(exp(gammaln(kl) -gammaln(m) -gammaln(kl-m+l)+m*log(p)+(kl-m) *log(l-p) ) ) ; 


777.77777.77.77.7.7777.77777.77777.777 


function  POv=mofn_steadystate (pOf ,h) 

7  POv=mofn_steadystate (pOf ,h) 

7.  Returns  the  steady-state  probabilities  of  each  state 
7  given  a  long-term  probability  of  succes  of  pOf  (e.g., 
7  false  alarm  probability)  and  assuming  no  successes. 

7. 


a 
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7.  Note  that  when  pOf  is  small,  POv  is  approximately 
/  POv  ~  [10  ...  0] 

7. 

if  p0f>0, 

P0v=p0f . ~binnum(h) ; 

P0v(end)=0;  70  Set  accepting  state  probability  to  zero 
P0v=P0v/sum(P0v) ;  7.  Condition  on  being  in  non-accepting  state 
else 

P0v=[l;  zeros(length(h)-l,l)] ; 
end; 


B-4  Subroutines  for  m-of-rt  detector  design 


function  [mv,nv, loglOTv] =mofn_f ar_spec(p0 , logl0T,nsol) 
7,  [mv , nv ,  loglOTv]  =mof n_f ar_spec  (pO ,  loglOT , nsol) 

7.  Returns  nsol  pairs  of  (m,n)  meeting  the  FAR  spec 
7.  and  the  value  of  loglO(Tbar)  achieved  by  the  pair 
icmax=20 ; 
if  nargin<3, 
nsol=5 ; 
end; 

m0=ceil (logl0T/-logl0 (pO) ) ; 

Iogl0e=logl0(exp(l)) ; 
for  i=l:nsol. 


m=m0+i-l ; 
mv ( i ) =m ; 
n0=mv(i) ; 


err=2 ; 
icv=0 ; 

while  (err>0.6)&&(icv<icmax)&&(n0<1000) , 


icv=icv+l ; 

f a=logl0e* (-gammaln(nO)+gammaln(nO-m+l)+gammaln(m) -m*log(p0) - (nO-m+1) *log(l-p0) ) ; 
fpa=logl0e* (-psi (nO) +psi (nO-m+1) -log ( 1-pO) ) ; 
deln=(fa-loglOT)/fpa; 
n0=n0-deln; 
err=abs (deln) ; 
end; 

nv(i)=f loor (nO) ; 
end; 

logl0Tv=logl0e* (gammaln(mv)+gammaln(nv-mv+l)-gammaln(nv) -mv*log(p0)- (nv-mv+1) *log(l-p0) ) ; 


77.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7. 
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function  [B,pd,k] =mofn_pd_breakpoint (pO,pl ,loglOT,Lmax) 

•/.  [B ,  pd ,  k]  =mof  n_pd_breakpoint  (pO ,  p  1 ,  loglOT ,  Lmax) 

•/.  Returns  the  breakpoints  for  m-of-n  detector  design  to  maximize 
•/.  the  lower  bound  on  Pd  formed  from  the  m-of-k  Pd  envelopes 

7. 

7,  pO  =  noise-only  success  probability 
7.  pi  =  noise+signal  success  probability 
7.  loglOT  =  loglO(Tbar)  FAR  specification 
7.  Lmax  =  maximum  signal  length  of  interest 

7. 

7.  B.mv  =  value  of  m  in  each  region 

7.  B.nv  =  maximum  value  of  n  meeting  FAR  spec  for  each  m 
7.  B.LO  =  beginning  of  each  region 
7.  B.L1  =  end  of  each  region 
7.  pd  =  Pd  lower  bound 

7.  k  =  vector  of  signal  lengths  associated  with  pd 

7.  Get  (m,n)  pairs  meeting  FAR  spec 
nsol=5 ; 

[mv , nv , loglOTv] =mof n_f ar_spec (pO , loglOT , nsol) ; 
while  max(nv) <Lmax+l , 
nsol=nsol+l ; 

[mv , nv , loglOTv] =mof n_f ar_spec (pO , loglOT , nsol) ; 
end; 

nsol=f ind (nv>Lmax , 1 , ’ f ir st ’ ) ; 
mv=mv(l :nsol) ; 
nv=nv(l :nsol) ; 
pd=zeros(l ,Lmax) ; 

[pdO ,k0] =mofn_pd_k_less_n(pl ,mv(l) ,nv(l)) ; 

B.mv(l)=mv(l)  ;  7»  currently  active  value  of  m 
B . nv ( 1 ) =nv ( 1 ) ; 

B.L0(1)=1;  7.  first  signal  length  for  this  region 

B.Llc(l)=nv(l) ; 

ir=l; 

pd(kO)=pdO; 

k_new=  []  ; 

for  i=l : (nsol-1) , 

[pdl ,kl] =mofn_pd_k_less_n(pl ,mv(i+l) ,nv(i+l)) ; 
if  pdl (end) >pdO (end) , 
ir=ir+l ; 

B.mv(ir)=mv(i+1) ; 

B.nv(ir)=nv(i+1) ; 

B . L0(ir)=f ind(pdl>pdO(end) , 1 , 1  first ’ ) ; 

B.Llc(ir)=nv(i+l) ; 

B.Ll(ir-l)=B.LO(ir)-l; 


UNCLASSIFIED 


approved  for  public  release 


56 


CausaSci  LLC 


UNCLASSIFIED 


2011-01 


%  Interim  region  is  flat  at  previous  Pd  (lower  bound) 
k_int=B . Lie (ir-1) :min( [Lmax  B.Ll(ir-l)] ) ; 
pd(k_int)=pdO (end) *ones(size(k_int) ) ; 

7«  Fill  in  up  to  end  of  new  value  of  m 
k_new=B.LO(ir) :min(Lmax,nv(i+l) ) ; 
pd(k_new)=pdl (k_new) ; 
pdO=pdl ; 
end; 
end; 

npd=length(pd) ; 
if  any(k_new), 

pd(k_new(end) :npd)=pdl (k_new(end) :npd) ; 
end; 

B.Ll(ir)=Lmax; 
k=l : length (pd) ;% 


mmxmmnmmramra 


function  [mO ,n0 , TO ,D0 , 10] =mof n_opt_mn_T_and_D (pO ,pl , loglOT) 
7.  [mO , nO ,  TO ,  DO ,  10]  =mof  n_opt_mn_T_and_D  (pO , pi ,  loglOT) 

7.  pO  =  noise-only  success  probability 
7,  pi  =  noise+signal  success  probability 

7.  loglOT  =  loglO(Tbar)  FAR  specification  (can  be  a  vector) 

7. 


7.  mO  =  optimal  value  of  m 

7.  nO  =  optimal  value  of  n 

7.  TO  =  achieved  average  time  between  false  alarms 
7.  DO  =  achieved  latency 

7.  10  =  structure  containing  initialization 


7,  Initial  starting  point 

[msO ,ns0] =mofn_opt_mn_T_and_D_init (pO ,pl , loglOT) ; 

10 .ms=ms0 ; 

I0.ns=ns0; 

for  i=l : length(loglOT) , 
fprintf ( ’ 

nsol=ms0(i) -ceil (loglOT (i) /loglO(pO) )+5 ; 

[ms ,ns , loglOTs] =mofn_f ar_spec (pO, loglOT (i) ,nsol) ; 
Dmin=inf ; 

isol=f ind(ms==msO(i) ) ; 

Dsl=mofn_stats_cf cn(pl ,ms(isol) ,ns (isol) ) ; 
I0.D(i)=Dsl; 

10 . T(i)=10~logl0Ts (isol) ; 

7o  Back  up  until  increasing 
Ds0=Dsl-l ; 
i0=isol-l ; 
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Ds2=  []  ; 

while  (iO>=l)&&(DsO<Dsl) , 

DsO=mofn_stats_cf cn(pl ,ms (iO) ,ns(iO)) ; 

7«  Initialization  beyond  minimum,  move  down  again 
if  DsCKDsl, 
isol=iO ; 
iO=iO-l; 

Ds2=Dsl ; 

Dsl=DsO ; 
end; 
end; 

•/.  Increase  until  at  minimum  if  haven’t  already  found  it 
if  ~any(Ds2) , 

Ds2=Dsl-l; 

while  (Ds2<Dsl)&&(isol<nsol) , 
l  Evaluate  D 

Ds2=mofn_stats_cf cn(pl ,ms (isol+1) ,ns(isol+l) ) ; 

7,  Not  at  minimum  yet,  keep  going 
if  Ds2<Dsl, 
isol=isol+l ; 

Dsl=Ds2; 

end; 

end; 

end; 

mO(i)=ms(isol) ; 
nO(i)=ns(isol) ; 

D0(i)=Dsl ; 

T0(i)=10. ~logl0Ts (isol) ; 
end; 


function  [ms ,ns] =mofn_opt_mn_T_and_D_init (pO,pl , loglOT) 
•/.  [ms  ,ns]  =mofn_opt_mn_T_and_D_init  (pO  ,pl ,  loglOT) 

7„  Returns  the  starting  point  design  (ms, ns)  to  minimize 
7„  Dbar  given  pO,  pi,  and  loglOT 


•/.  Get  efficiency  from  Page’s  test  &  Dbar 
bl=l/(l+log(pl/pO)/log((l-pO)/(l-pl))); 

[“ , ~ , eta] =page_stats_bino (pO ,pl ,bl , 1) ; 

Dp=loglOT/(eta*loglO(exp(l))) ; 
for  i=l : length(loglOT) , 

%  Get  Pd  breakpoints 

[B , ~ , ~] =mof n_pd_breakpoint (pO ,pl , loglOT (i) , 5*ceil (Dp (i) ) ) ; 

7,  Find  m  maximizing  Pd  for  signal  of  length  Dbar  from  Page’s  test 
iopt=f ind(B .  L0<=Dp(i) , 1 , ’ last ’ ) ; 
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ms(i)=B.mv(iopt) ; 
ns(i)=B.nv(iopt) ; 
end; 


mnmnmmy.nmnmm 


function  [T,D,eta] =page_stats_bino (pO ,pl ,b,h) 

7„  [T,D, eta]  =page_stats_bino(pO,pl  ,b,h) 

•/.  Returns  performance  measures  for  Page’s  test 

7„  operating  on  binary  data 

“/«  pO  =  noise-only  success  probability 

“/«  pi  =  noise+signal  success  probability 

•/.  b  =  bias 

V.  h  =  threshold 

T0L=le-4; 

°/0  Get  unity  roots  of  MGF 
tO=mgf_root_bino(pO,b,TOL) ; 
tl=mgf _root_bino (pi ,b,T0L) ; 

°/0  Approximations 

ugO=pO-b; 

ugl=pl-b; 

T=nan*ones (size (h) ) ; 

D=T ; 

ib=f ind(h>l-b) ; 

T(ib)=(l+h(ib)*tO-exp(h(ib)*tO)) ./(tO*ugO) ; 
D(ib)=(l+h(ib)*tl-exp(h(ib)*tl)) ./(tl*ugl) ; 
eta=tO*ugl ; 


function  t=mgf _root_bino (p , b , TOL) 

“/«  t=mgf_root_bino(p,b,TOL) 

•/.  Returns  the  unity  root  of  the  moment  generating  function  (MGF) 

“/«  for  binary  data  with  success  probability  p  and  bias  b 
if  nargin<3, 

T0L=le-6 ; 
end; 

if  (p-b)<0, 

t=-log(p)/ (1-b) ; 
else 

t=log(l-p)/b; 

end; 

err=2*T0L ; 
while  err>T0L, 

f =log(p*exp( (1-b) *t)+(l-p) *exp(-b*t) ) ; 

dfdt=(p*(l-b)*exp(t*(l-b))-b*(l-p)*exp(-t*b))/ (p*exp((l-b)*t)+(l-p)*exp(-b*t)) ; 
t=t-f/dfdt ; 


UNCLASSIFIED 


approved  for  public  release 


59 


CausaSci  LLC 


UNCLASSIFIED 


2011-01 


err=abs(f/dfdt) ; 
end; 
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